Doctorat De L'universite

D E Lyon

Devant Le

M Laquerriere

Rapporteur Mme Vincent

M Douek

M Desvignes Michel

Ecole Sigle

Doctorale

Renée El

Melhem Bât

Pascal Blaise

R Gourdon

M Stéphane

ÉLECTROTECHNIQUE E E A Électronique

M Gérard Scorletti
email: gerard.scorletti@ec-lyon.fr

M Philippe
email: philippe.normand@univ-lyon1.fr

Mme Emmanuelle
email: emmanuelle.canet@univ-lyon1.fr

Étage

M Luca

Jean-Yves Buffière
email: jean-yves.buffiere@insa-lyon.fr

M Jocelyn Bonjour

Viviane Polsinelli
email: viviane.polsinelli@univ-lyon2.fr

Brigitte Dubois

J Y Toussaint

M Christian Montes

N°d'ordre NNT : 2018LYSEI056

Keywords: Computed tomography, X-ray, Spectral CT, Photon counting detector, material decomposition, simulation Tomographie, Rayons X, CT spectrale, Détecteur à comptage des photons, Décomposition de matériaux, Simulation

X-ray computed tomography (X-ray CT) plays an important part in non-invasive imaging since its introduction. During the past few years, numerous technological advances in Xray CT have been observed, including spectral CT, which uses photon counting detectors (PCDs) to discriminate transmitted photons corresponding to selected energy bins in order to obtain spectral information with one single acquisition.

Spectral CT enables us to overcome many limitations of the conventional CT techniques and opens up many new application possibilities, among which quantitative material decomposition is the hottest topic. A number of material decomposition methods have been reported and dierent experimental systems are under development for spectral CT. According to the type of data on which the decomposition step operates, we have projection domain method (decomposition before reconstruction) and image domain method (decomposition after reconstruction).

The commonly used decomposition is based on least square criterion, named proj-LS and ima-LS method. However, the inverse problem of material decomposition is usually illposed and the X-ray spectral CT measurements suer from Poisson photon counting noise.

The standard LS criterion can lead to overtting to the noisy measurement data. In the present work, we have proposed a least log-squares criterion for projection domain method to minimize the errors on linear attenuation coecient: proj-LLS method. Furthermore, to reduce the eect of noise and enforce smoothness, we have proposed to add a patchwise regularization term to penalize the sum of the square variations within each patch for both projection domain and image domain decomposition, named proj-PR-LLS and ima-PR-LS method.

The performances of the dierent methods were evaluated by spectral CT simulation studies with specic phantoms for dierent applications: (1) Medical application: iodine and calcium identication. The decomposition results of the proposed methods show that calcium and iodine can be well separated and quantied from soft tissues. (2) Industrial application: ABS-ame retardants (FR) plastic sorting. Results show that 3 kinds of ABS materials with dierent ame retardants can be separated when the sample thickness is favorable.

Meanwhile, we simulated spectral CT imaging with a PMMA phantom lled with Fe, Ca and K solutions. Dierent acquisition parameters, i.e. exposure factor and number of energy bins were simulated to investigate their inuence on the performance of the proposed methods for iron determination.

I am grateful for the past three and a half years of my Ph.D. study in the laboratory LVA at INSA de Lyon. I have benefited a lot from the excellent research environment in the laboratory. I deeply believe that this thesis will not have been finished without the support and kindness of my supervisors, my colleagues, my friends and my family.

I would like to express my sincere gratitude to my supervisors Valérie Kaftandjian, Philippe Duvauchelle and Yuemin Zhu. Their insightful and inspiring guidance, their continuous support and encouragement are indispensable to the accomplishment of this thesis. I learned a lot from their serious attitude and passion for life and work. Great thanks should be given to Valérie who has devoted great time, patience and efforts on the guidance of my work.

I also would like thank Prof. Patrice Laquerriere and Nicole Vincent for reviewing my thesis and Prof. Philippe Douek, Michel Desvignes and Su Ruan for accepting to be the members of the committee. Thank you for all the comments and questions that help to improve the quality of this thesis.

I would like to thank prof. Jérome Antoni who provided insight and expertise through fruitful discussions. I thank all my colleagues in the laboratory for their helps and friendliness during my life in France. I will always remember the precious time we have spent together. I must thank my close and dear friends who keep me company and gave me the courage on the road of PhD study.

I sincerely thank my grandparents, my father, my sisters and all my family for their endless love and support in my life. I want to especially thank my boyfriend Shishun for being very understanding and considerable, for the beautiful memories of travelling and for the countless hours we have spent on the TGV between Rennes and Lyon every weekend.

Finally, I would like to dedicate this thesis to my dear mother in heaven. You were a big part of my life and made me who I am today. I will miss you and love you as always.

Introduction

La tomographie (CT -computed tomography) aux rayons X est introduite pour la première fois dans les années 1970 par Cormack et Hounseld, pour laquelle ils ont reçu le prix Nobel de physiologie ou de médecine en 1979. Cette invention a montré la possibilité de donner des images en coupe de haute qualité avec des mesures multi-angles et des calculs mathématiques complexes. La Figure 1 montre la mécanique de l'imagerie CT, dont les composants de base sont la source de rayons X, l'objet et le détecteur. Par la suite, l'avancement de la CT aux rayons X a grandement bénécié des algorithmes de reconstruction développés, de la technologie informatique avancée et des améliorations technologiques des sources et détecteurs. Il y a maintenant de nombreuses applications dans les domaines médical, industriel, caractérisation des matériaux, aérospatial et géologique. Le concept de CT résolue en énergie a été discuté depuis l'invention de la CT. Dans la publication référence [Hounseld, 1973] de Hounseld, il a déjà été proposé de prendre des mesures de la même tranche d'un matériau à deux tensions diérentes an de déterminer le numéro atomique des matériaux. En 1976, Alvarez et Macovski ont proposé d'extraire l'information dépendant de l'énergie en représentant l'atténuation totale avec des fonctions correspondant à l'eet photoélectrique et Compton, et les parties indépendantes de l'énergie des deux eets ont ensuite été reconstruites respectivement [START_REF] Alvarez | Energy-selective reconstructions in x-ray computerised tomography[END_REF]. En 1977, Riederer et al. ont présenté leur travail d'imagerie sélective de l'iode par CT en utilisant trois faisceaux de rayons X fortement ltrés [START_REF] Riederer | Selective iodine imaging using k-edge energies in computerized x-ray tomography[END_REF]. L'imagerie sélective des matériaux utilisant la CT à bi-énergie (DECT dual energy computed tomography) s'étend également à l'imagerie radiographique aux rayons X [START_REF] Brody | Dual-energy projection radiography: initial clinical experience[END_REF], Lehmann et al., 1981] désormais couramment utilisée dans les systèmes de contrôle de bagages dans les aéroports.

Suivant le concept ci-dessus, une variété de systèmes DECT a été construite pour l'imagerie à résolution d'énergie. Il existe actuellement trois techniques principales pour réaliser la DECT. Le premier est appelé technique de commutation rapide des kilovolts, où la tension du tube à rayons X est rapidement commutée pour produire deux faisceaux de rayons X avec des énergies diérentes pour balayer toutes les autres projections pendant toute la durée du balayage. Son principe est représenté sur la Figure 2 (b). Le second est une technique CT à double source dans laquelle deux tubes à rayons X et deux détecteurs correspondants sont équipés dans le système pour acquérir simultanément deux ensembles de données pour la même tranche de l'objet. Ces deux techniques utilisent toutes deux des spectres de rayons X diérents et leurs principes sont représentés sur la Figure 2 (b). La troisième est une technique de détection à double couche (détecteur sandwich), qui utilise encore un détecteur intégrant l'énergie, mais deux ensembles de réseaux de scintillateurs et de photodiodes empilés les uns sur les autres, de sorte que les photons de basse énergie sont absorbés par le réseau scintillateur supérieur, tandis que les photons de plus haute énergie traversent la couche supérieure et sontabsorbés par le réseau de scintillateurs du dessous [Nasirudin, 2015], voir la Figure 2 (c). A la diérence de la DECT, la CT spectrale emploie un détecteur à comptage de photons (PCD photon counting detector) qui est capable d'obtenir des informations spectrales de plusieurs bandes d'énergie avec une seule acquisition, Figure 2 (d). La CT spectrale permet de surmonter de nombreuses limitations des techniques précédentes et ouvre de nombreuses nouvelles applications. La décomposition quantitative des matériaux est le sujet le plus étudié, y compris l'imagerie K-edge pour des matériaux de contraste et l'imagerie sélective d'autres matériaux auparavant indiscernables. D'autres applications incluent la pondération énergétique et l'imagerie monochromatique virtuelle pour améliorer la qualité de l'image.

Dans le présent travail, nous nous concentrons sur l'application de la décomposition quantitative des matériaux par CT spectrale. Un certain nombre de méthodes de décomposition des matériaux ont été rapportées et diérents systèmes expérimentaux sont en cours de développement pour la CT spectrale. Cependant, la plupart des recherches se concentrent sur l'imagerie des agents de contraste avec des numéros atomiques élevés (Z ), tels que l'iode, le gadolinium et l'or, qui ont des raies K uniques dans la gamme d'énergie de détection. Il est très dicile de distinguer les éléments qui ont de faibles numéros atomiques, surtout quand leurs nombres Z sont proches. Dans cette phase primaire des systèmes ex- périmentaux, il est dicile d'obtenir les données brutes pour le traitement de la recherche.

La simulation est un moyen ecace, rentable et réalisable de résoudre ce problème. Par conséquent, la thèse se concentre autour de trois aspects: (1) l'étude des méthodes de décomposition en matériaux pour la CT spectrale; (2) 2 Méthodes de décomposition en matériaux

Selon la formulation de la décomposition en matériaux, les méthodes de décomposition existantes peuvent être divisées en deux catégories. La première est la décomposition basée sur l'eet physique, également appelée imagerie K-edge, où l'atténuation est considérée comme une combinaison linéaire de l'eet photoélectrique, de l'eet Compton et des composants K-edge s'il y a des matériaux dont le saut d'absorption K est compris dans la gamme d'énergie du détecteur de CT. La deuxième est la méthode basée sur les matériaux, notre travail était autour de celle-ci, où le coecient d'atténuation linéique est décrit comme les contributions d'une base de matériaux. Ensuite, selon le type de données sur lequel l'étape de décomposition fonctionne, les méthodes peuvent également être divisées en méthodes fonctionnant dans le domaine des projections (sinogramme ou projection radiographique) ou bien dans le domaine de l'image reconstruite (coupe 2D ou volume 3D). Le premier type, également appelé méthode de pré-reconstruction, décompose les données de projection en plusieurs composants désirés (i.e. base des matériaux), puis applique une reconstruction indépendante à chaque composant pour obtenir leur distribution respective dans l'espace. Le deuxième type, également connu sous le nom de méthode de post-reconstruction, fonctionne sur des images CT reconstruites de manière conventionnelle correspondant aux diérentes bandes d'énergie acquises pour obtenir la concentration de matériaux cibles en utilisant une inversion de matrice linéaire, voir la Figure 3. Les deux méthodes ont leurs avantages et leurs inconvénients, en fonction de diérentes situations.

Dans la suite, ces deux méthodes seront dénommées projection ou image et repérées par le préxe proj ou ima. 3 Résultats de décomposition en matériaux pour diérentes applications: étude de fantômes L'objectif de cette partie du travail est d'introduire la méthode de simulation de l'imagerie par rayons X et d'évaluer les méthodes de décomposition dans le domaine des projections proposées à travers l'étude de simulation. Nous simulons le système CT / CR spectral qui possède des capacités de résolution en plusieurs bandes énergétiques en utilisant le logiciel Virtual X-ray Imaging (VXI) [START_REF] Duvauchelle | A computer code to simulate x-ray imaging techniques[END_REF]. Ce logiciel a été développé à l'INSA de Lyon (France) pour simuler l'imagerie radiographique, radioscopique et tomographique. La simulation est basée sur des techniques de tracé de rayons et s'appuie sur la loi d'atténuation des rayons X de Beer-Lambert. Ciblant une application médi-cale (imagerie de l'athérosclérose, CT) et une application industrielle (tri des plastiques, radiographie (CR)), nous simulons à travers le logiciel VXI le processus d'imagerie par rayons X avec les deux fantômes correspondants, nous appliquons les méthodes proposées et analysons leurs performances. Pour la première application, nous allons comparer les performances de décomposition de la méthode proj-LS, proj-LLS et proj-PR-LLS et pour la deuxième application, nous utiliserons la méthode proj-PR-LLS que nous complétons par une étape de classication (dénommée proj-PR-LLS-clas). Pour les deux applications de ce chapitre, la CT spectrale ou le processus d'imagerie CR sont simulés avec les paramètres d'acquisition habituellement utilisés dans l'application considérée, par ex. nombre de bandes d'énergie et facteur d'exposition.

Nous comparerons la performance de la méthode proj-PR-LLS et ima-PR-LS, ainsi que l'inuence des paramètres d'acquisition sur leur performance dans le prochain chapitre.

3.1 Application médicale: identication de l'iode et du calcium

Les maladies cardiovasculaires (MCV) engendrent d'énormes conséquences sanitaires et économiques à l'échelle mondiale et l'athérosclérose est un contributeur majeur des MCV.

L'athérosclérose est une maladie dans laquelle l'artère interne se rétrécit en raison de dépôt de plaque. La présence de plaques à l'intérieur des artères coronaires est un facteur de risque grave pour les événements cardiovasculaires indésirables. Les ruptures de plaques vulnérables sont à l'origine d'environ 70% des infarctus aigus du myocarde et des décès coronariens soudains. Il est très important de détecter la présence de plaque et d'évaluer sa vulnérabilité.

Les méthodes existantes peuvent détecter la plaque à haut risque dans une certaine mesure avec des capacités limitées et une évaluation améliorée peut être réalisée par une combinaison de plusieurs méthodes. Cependant, il serait avantageux de trouver une modalité d'imagerie capable d'évaluer en même temps plusieurs caractéristiques de la plaque à haut risque et de permettre une identication avec une meilleure précision.

Pour l'inspection de la plaque, l'identication de ses composants est d'une importance signicative. Sachant que la CT spectrale aux rayons X possède une capacité de décomposition quantitative en matériaux comparée à l'angiographie conventionnelle par CT, nous étudions la capacité de décomposition des méthodes proposées (méthode proj-LS, proj-LLS et proj-PR-LLS) et évaluons leur performance dans l'identication de l'iode, du calcium et des tissus mous au moyen d'une étude sur le fantôme thoracique humain.

Fantôme simulé & géométrie du système

Nous utilisons un fantôme thoracique humain pour simuler les situations in vivo. Les composants et les densités de diérents tissus inclus dans le fantôme du thorax (Figure 5 (a)) sont dénis selon les publications de White [START_REF] White | Average softtissue and bone models for use in radiation dosimetry[END_REF] et de la publication 44 de la Commission internationale sur les unités et mesures radiologiques (ICRU44) [ICRU, 1989]. En plus des diérents types de tissus corporels, des inserts supplémentaires de

Résultats de la décomposition

Dans cette application, nous avons choisi le tissu mou, le calcium et l'iode comme base de matériaux pendant la décomposition. Trois méthodes de décomposition de matériaux dans le domaine des projections avec diérentes fonctions de coût ont été appliquées: proj-LS, proj-LLS et proj-PR-LLS. La gure 7 montre les images décomposées en utilisant les trois méthodes. On observe que toutes les méthodes donnent la distribution des tissus mous, du calcium et de l'iode:

• La peau, les tissus adipeux, les tissus mous moyens, le cartilage, l'aorte, les poumons et le c÷ur ne sont visibles que dans l'image des tissus mous. En outre, une partie des os (côtes, vertèbres et sternum) est également présente dans l'image des tissus mous, car les composants des os contiennent non seulement du calcium, mais aussi d'autres éléments comme l'hydrogène et l'oxygène.

• L'image du calcium met bien en évidence les os et les solutions de calcium insérées.

• L'iode est séparée de manière appropriée dans l'image spécique à l'iode.

Page x

Cette thèse est accessible à l'adresse : http://theses.insa-lyon.fr/publication/2018LYSEI056/these.pdf © [T. Su], [2018], INSA Lyon, tous droits réservés

Conclusion and discussion

Motivés par la nécessité de distinguer les compositions de la plaque athérosclérotique, nous avons simulé le processus d'imagerie par CT spectrale avec un fantôme thoracique humain réaliste, et étudié les performances de décomposition en matériaux des trois méthodes présentées dans le chapitre précédent: proj-LS, proj-LLS et proj -PR-LLS.

ABS basis

With noise; r=10 4.4 . 6), (0.1, 10), (0.1, 15), (0.1, 20), (0.1, 30), (0.1, 60), puis [START_REF]Problem statement and objective[END_REF][START_REF][END_REF], [START_REF]Problem statement and objective[END_REF]10), [START_REF]Problem statement and objective[END_REF]15), et ainsi de suite. L'axe Y représente diérentes concentrations de fer. Les résultats de décomposition avec un facteur d'exposition de 0.025 mA•s ne sont pas présentés dans ce graphique pour des raisons de commodité de lecture.

2.2

The electromagnetic spectrum, presented as a function of wavelength, frequency, and energy, retrieved from [Seibert, 2004]. . . . . . . . . . . . . . .

2.3

During the photoelectric event, the incident photon is annihilated. Part of its energy is consumed by ionization of atom and the rest is given to ejected photoelectron as its kinetic energy. The transfer of electron from outer shell to inner shell to ll in the vacancy created by photoelectron will generate either a. X-ray uorescence, or b. Auger electron. Figure is inspired by [START_REF] Barrett | Radiological imaging: the theory of image formation, detection, and processing[END_REF]. . . . . . . . . . . . . . . . . . . . .

2.4

The mechanism of Compton scattering: incident photon is scattered, its energy and direction change and a recoiling electron is ejected. Figure is inspired by [START_REF] Barrett | Radiological imaging: the theory of image formation, detection, and processing[END_REF]. . . . . . . . . . . . . . . . . . . . .

2.5

The mechanism of pair production: incident photon is annihilated, an electronpositron pair is generated. Figure is inspired by [START_REF] Barrett | Radiological imaging: the theory of image formation, detection, and processing[END_REF]].

2.6

Energy dependence of linear attenuation coecient (water): a. photoelectric absorption µ pe ; b. Compton scattering µ C ; c. pair production µ pp ; d. total attenuation coecient µ tot . Points e1 and e2 are energies where µ pe = µ C and µ C = µ pp . Data taken from database XGAM [START_REF] Berger | Nist x-ray and gamma-ray attenuation coecients and cross sections database[END_REF].

2.8

The mechanics of CT imaging: cross sectional image is reconstructed from projection data measured from multiple angles. . . . . . . . . . . . . . . . .

2.9

The EMI Mark I scanner (a), and the cross sectional image of brain (b). . .

3.5

The material decomposition methods discussed and proposed (in blue) in this chapter. The term CR in refers to computed radiography. . . . . . .

4.1

Illustration of the simulation principles. The ray SK intersects two meshes at points A and B. Geometrical calculations enable the attenuation path length AB to be determined. Ray (1): transmitted photons. Ray ( 2) and

(3): scattered photons. Figure is retrieved from [START_REF] Duvauchelle | A computer code to simulate x-ray imaging techniques[END_REF]. .

5.3

Mass attenuation coecients (µ m ) of the 5 materials (Fe, Ca, K, PMMA and water) contained in this phantom within the detecting energy range from 30 keV to 90 keV. Data taken from [START_REF] Hubbell | Tables of x-ray mass attenuation coecients and mass energy-absorption coecients (version 1.4)[END_REF]. . . . . List of Tables 2.1 Comparison of dierent PCDs. This table is not a complete list of all available PCDs and is only for comprehensive review. . . . . . . . . . . . . . . . It is known that the conventional CT with single X-ray tube and energy integrating detector measures the total attenuation of an object over the entire spectrum. However, various compositions of object respond dierently to a particular energy of the incoming radiation. Making use of such property, we can distinguish dierent components of the object. In the past years, numerous technologies of X-ray CT have been proposed, including dual-energy CT (DECT) and spectral CT. DECT has limited ability in material decomposition by producing two data sets corresponding to two dierent x-ray spectra. In contrast, spectral CT employs photon counting detectors (PCDs) to obtain more complete spectral information: the transmitted photons can be discriminated corresponding to selected energy bins with one single acquisition. In this way, multiple materials can be simultaneously identied. A number of material decomposition methods have been reported and dierent experimental systems are under development for spectral CT. However, most researches focus on the imaging of contrast agents with high atomic numbers (Z number), such as iodine, gadolinium and gold, which have unique K-edges within detection energy range.

It is much more dicult to distinguish those having low atomic numbers, especially when their Z numbers are close.

Some of the material decomposition methods proposed for spectral CT imaging can also be applied to PCD-based spectral computed radiography (CR). Targeting on one medical application of spectral CT (atherosclerosis imaging) and one industrial application of spectral CR (plastic sorting), our objective is to develop material decomposition methods to improve the imaging quality of these two applications:

• Medical application: atherosclerosis imaging. Atherosclerosis is a disease in which the inside artery narrows due to the build-up of plaques. Plaques inside the coronary arteries are a serious risk factor for adverse cardiovascular events. The ruptures of vulnerable plaques are the cause of about 70 % of fatal acute myocardial infarctions and sudden coronary deaths. Discrimination of plaque compositions is very important in the identication of high-risk plaques.

• Industrial application: plastic sorting. Waste electrical and electronic equipments (WEEE) have been increasing rapidly due to the development of electronic industry. They contain various compositions among which polymers occupy an important part. Due to the presence of electronic power of electrical and electronic equipment (WEEE), the plastic materials should meet high re safety standards. However, it is not possible to realize with pure polymers, therefore ame retardant (FR) are added to change ammability of plastics and increase the re resistance. Recycling of plastics from WEEE is challenging because of the existence of FRs (especially the brominated and chlorinated FRs), which might result in serious environmental pollution. Therefore, the sorting technique is important for subsequent treatment and recycling of WEEE to avoid major environmental and health problems.

At this primary stage of the proposed new research, it is hard to obtain the experimental data due to the non-existence of physical or commercial spectral CT systems. Simulation is then an ecient, cost-eective and feasible way to address the problem. Therefore, we will study the material decomposition methods, evaluate their performance through simulations in view of the above-mentioned applications, and further investigate the inuence of dierent acquisition parameters on material decomposition performance.

Main contributions

Performance evaluation The main contributions of this thesis include: the proposal of new material decomposition methods, the validation and assessment of the proposed methods through simulated spectral CT or CR imaging, and the investigation of the inuence of dierent acquisition parameters on material decomposition performance, as illustrated in Figure 1.1. In the following, we detail the contributions following the three aspects.

Plastic sorting

Atherosclerosis imaging

Proposed material decomposition methods

Projection domain methods

Image domain methods

Investigation of the influence of different acquisition parameters on material decomposition performance

• Material decomposition methods (Chapter 3).

-General understanding. Making use of multi-energy information, spectral CT has the potential to discriminate dierent components inside the object. Many material decomposition methods have been proposed and dierent experimental spectral CT systems are under development. According to the type of data on which the decomposition step operates, the methods can be divided into projection domain method and image domain method. The former, also called pre-reconstruction decomposition method, rstly decomposes projection data 1.2. MAIN CONTRIBUTIONS into several desired components (i.e. materials), and then applies independently reconstruction to each component to obtain their respective distribution in space. The latter, also known as post-reconstruction method, operates on conventionally reconstructed CT images corresponding to dierent energy bins and it obtains the concentration of target materials using linear matrix inversion.

-Proposed methods. The commonly used objective function to solve this problem is based on least squares criterion, which minimizes the errors between expected and measured number of photons (projection domain method) or the errors between expected and reconstructed attenuation coecients (image domain method). However, the inverse problem of material decomposition is usually ill-posed and the X-ray spectral CT measurements suer from Poisson photon counting noises. The standard least squares objective function can lead to overtting to the noisy measurement data. Therefore, we have proposed a least log-squares criterion for projection domain method to minimize the errors on linear attenuation coecient. Furthermore, we have proposed to add a regularization term to penalize the sum of the square variations of the decomposed images, in order to reduce the eect of noise and enforce smoothness. The proposed patchwise regularized method is applied to both projection domain and image domain decomposition.

• Validation and assessment via simulation (Chapter 4).

-Simulation of spectral X-ray imaging. Virtual X-ray imaging (VXI) software is used to simulate the spectral X-ray CT and CR imaging of specic phantoms.

The simulation is based on ray-tracing techniques together with the X-ray attenuation law.

-Medical application. A realistic computational human thorax phantom has been built for the simulation of spectral CT imaging. Three projection domain methods, including the proposed ones, have been applied for the identication of soft tissue, iodine and calcium. Their performance are compared in terms of the average error and contrast-to-noise ratio.

-Industrial application. The spectral CR imaging of a phantom with multiple cubes containing dierent plastic materials has been simulated. We have used the proposed decomposition methods to identify the type of ame retardant contained in each cube in order to achieve plastic sorting.

• Inuence of dierent acquisition parameters (Chapter 5).

We have investigated the inuence of dierent acquisition parameters, i.e., exposure factor and the number of energy bins, on material decomposition performance through a phantom study for iron determination. We have simulated spectral CT imaging of the phantom with dierent settings of exposure factor (0.025, 0.1, 1, 10 and 100 mA•s per projection) and number of energy bins (6, 10, 15, 20, 30 and 60),

and have compared the material decomposition performance of image domain and projection domain methods with the simulated data.

The thesis manuscript is organized as follows:

In Chapter 2, entitled General principles of spectral CT, the principles of the interactions of X-rays with matter are introduced, including photoelectric absorption, Compton scattering and pair production; the principles of X-ray CT are described, including conventional energy integrated CT, dual-energy CT and spectral CT.

In Chapter 3, entitled Material decomposition methods, a review of the existing material decomposition methods for spectral CT is presented. Based on the existing methods, we have proposed patchwise regularized decomposition methods in both projection domain and image domain.

In Chapter 4, entitled Material decomposition results for dierent applications: phantoms study , the X-ray spectral imaging simulation procedure is detailed and the decomposition results of proposed methods are demonstrated through phantom study for both medical and industrial applications.

In Chapter 5, entitled Exploring the inuence of dierent acquisition parameters, a series of spectral CT acquisitions have been simulated with dierent parameters and the inuence of these parameters on material decomposition performance of the proposed methods is investigated.

In Chapter 6, entitled Conclusion and perspectives, a brief summary of the main results, the conclusions and future perspectives are presented.

In 1895, the German physicist Röntgen was investigating the glow that occurred during electric discharges inside an evacuated glass tube. During his experiments, he was surprised to notice that a screen coated with crystals of barium platinocyanide started to glow. The screen happened to be in his laboratory for the detection of ultraviolet radiation. In the following days, he carried out a series of tests and nally conrmed the existence of the mysterious rays [Cierniak, 2011]. He published his discovery in 1896 [Röntgen, 1896] with an enclosed X-ray picture of the hand of his wife Bertha (Figure 2.1), which clearly showed her bones in ngers and the wedding ring. His publication showed the possibility of noninvasive imaging of internal features and soon spread over the world. A lot of research was carried out on X-rays quickly and devices were built to take X-ray pictures. To reward this remarkable discovery, Röntgen was awarded the rst Nobel Prize in Physics. The electromagnetic spectrum, presented as a function of wavelength, frequency, and energy, retrieved from [Seibert, 2004].

The X-ray imaging approaches, including X-ray radiography and CT are based on the principle that the primary beam attenuates when it passes through the object. The number of photons decreases progressively in this process. According to Beer-Lambert law, the attenuation of a monochromatic X-ray beam in a homogeneous medium can be expressed as:

Φ(x) = Φ 0 exp(-µx), (2.1)
where the Φ 0 is the primary X-ray photon uence, Φ(x) is the attenuated photon uence with respect to traveling distance of x in the medium, µ is the linear attenuation coecient.

Another important parameter mass attenuation coecient µ m , is dened by:

µ m = µ/ρ, (2.2)
where ρ is material density. µ m has unit of centimeters squared per gram. Each element has its unique mass attenuation coecient curve, making it possible to distinguish dierent components of a mixture object.

There are three principal X-ray interactions with matters that account for the total attenuation. They are photoelectric absorption µ pe , Compton scatteringµ C , and pair production µ pp [ [START_REF] Barrett | Radiological imaging: the theory of image formation, detection, and processing[END_REF]. In fact, there is another eect that occurs when the X-ray interacts with the matter: Rayleigh scattering, also known as coherent scattering. It is an elastic interaction between a photon and an atom, where the photons direction is changed by a small angle but its energy is not. The scattered photons are mainly in the forward direction and the probability of this event occurring is low. Therefore it has little eect on the attenuation coecient. Thus the overall attenuation coecient µ tot and µ tot m can be considered as the sum of photoelectric absorption, Compton scattering and pair production:

µ tot = µ pe + µ C + µ pp , (2.3)
and

µ tot m = µ pe m + µ C m + µ pp m .
(2.4)

The following sections will introduce the detailed mechanism of photoelectric absorption, Compton scattering and pair production. [START_REF] Barrett | Radiological imaging: the theory of image formation, detection, and processing[END_REF].

In photoelectric eect, the incident photon with energy of hν 0 interacts with a target atom, its energy is totally absorbed by the atom by ionization and imparting kinetic energy E k to the ejected photoelectron:

hν 0 = E k + I, (2.5)
where I is the ionization potential, also called binding energy needed for a particular electron to be involved in a photoelectric event. Photoelectric interaction will occur only when the incident photon has energy that is equal to or greater than I. A photon having an energy just above the binding energy of the electron is more likely to be absorbed than a photon having an energy just below this binding energy, causing a sudden change in the attenuation coecient of the atom, called edge". A K-shell interaction is approximately four to ve times more probable than an L-shell interaction if both interactions are energetically allowed.

As analyzed in the above, when the photoelectric event occurs, a vacancy is created in the particular shell, which makes the electronic structure of the atom unstable. The electron from outer shell will transfer to ll this vacancy and create characteristic X-ray Page 10

Cette thèse est accessible à l'adresse : http://theses.insa-lyon.fr/publication/2018LYSEI056/these.pdf © [T. Su], [2018], INSA Lyon, tous droits réservés lines. The energy of the line is equal to the dierence in binding energies of the two shells. This will generate either X-ray uorescence or an Auger electron (by transferring the energy of X-ray line to another electron which is then emitted).

The linear attenuation contributed by photoelectric absorption µ pe is given by experiments as:

µ pe ≈ k Z m (hν 0 ) n ρ A , (2.6)
where k is a constant depending on the shell involved, Z is the atomic number and A is the atomic weight, m and n are parameters that varies slowly with Z and ν 0 . A useful estimation for parameter m and n is given as m = 4 and n = 3, making Equation (2.6) into

µ pe = k Z 4 (hν 0 ) 3 ρ A .
(2.7)

We can see that µ pe increases with atomic number and material density and is negatively related with energy of incident photon.

Compton scattering

Incident photon

Scattered photon

Recoiling electron

Original position of electron With α = hν 0 /m 0 c 2 (m 0 c 2 = 0.511meV), the scattered photon with the remaining energy hν ′ can be given by:

hν ′ = hν 0 (1 + α -αcosθ) -1 . (2.8) Page 11
Cette thèse est accessible à l'adresse : http://theses.insa-lyon.fr/publication/2018LYSEI056/these.pdf © [T. Su], [2018], INSA Lyon, tous droits réservés

INTERACTIONS OF X-RAYS WITH MATTER

Collision cross section σ is dened as the cross sectional area presented by the attenuating material to the incident beam. It describes the probability that an incident photon will experience a collision with material, and furthermore be attenuated. The relationship between collision cross section and linear attenuation coecient can be expressed as: σ = µ/n, n is the number of attenuating particles per cubic centimeter. µ = σn.

(2.9) Klein and Nishina rstly deduce the formula of collision cross section for Compton scattering σ C :

f KN = σ C σ 0 = 3 4 [ 2(1 + α) 2 α 2 (1 + 2α) + ln(1 + 2α) α ( 1 2 - 1 + α α 2 ) - 1 + 3α (1 + 2α) 2
] ,

(2.10)

where f KN is Klein and Nishina function, σ 0 = 8πr 2 0 is the cross section for classical Thomson scattering, r 0 = e 2 /(m 0 c 2 ) = 2.818 × 10 -13 cm.

The electron density n e (cm -3 ) in terms of Avogadro's number N 0 (6.02 × 10 23 mol -1 ) is:

n e = N 0 ρZ/A. (2.11)
Assume that all atomic electrons take part in the Compton scattering, then the linear attenuation coecient caused by Compton eect is:

µ C = n e σ C = (σ 0 N 0 ρZ/A)f KN .
(2.12)

It can be seen that Compton scattering correlates mainly with the density of materials, while the type of atom, especially when the atomic number is low, has little impact since Z/A ≈ 0.5. [START_REF] Barrett | Radiological imaging: the theory of image formation, detection, and processing[END_REF].

Pair production

In a pair production eect, the incident photon interacts with the nucleus and convert its energy into an electron-positron pair. Suppose the kinetic energy of the electron and positron as E e + and E e -, the photon must have higher energy than the sum of the rest mass energies of the electron and positron (2 × 0.511Mev = 1.022Mev), since:

hν 0 = E e + + E e -+ 1.022Mev.
(2.13)

According to the complicated calculations thorough quantum electrodynamics in the form of Feynman diagrams, the collision cross section for pair production σ pp is approximately positively related with Z 2 , and it increases with the incident photon energy.

Total attenuation coecient

From previous analysis, the total attenuation coecient can be considered sum of those three major interactions. It depends on the density, atomic number of materials, and also depends on the energy of incident photon. [START_REF] Berger | Nist x-ray and gamma-ray attenuation coecients and cross sections database[END_REF].

For each atom, there is a particular energy where µ pe = µ C and another energy where µ C = µ pp (Figure . 2.6, point e1 and point e2). By plotting the relationship of atomic number and those particular energies, we can get a general view of the dominant regions of each attenuating process, see Figure.

Conventional X-ray computed tomography 2.2.1 Historical background

After X-ray has been discovered, numerous research has been carried out for seeking new imaging modality. Many eorts have been done on planar X-ray radiography to produce 2D projections of objects situated between X-ray source and the detector. However, it is dicult to distinguish the overlapping components caused by the projection, such as the overlapping soft tissues and bones in physical scan. The introduction of X-ray computed tomography (CT) turned to be the right solution for this problem.

Based on the theory that materials have dierent attenuation properties to X-ray, CT scanners take measurements of integral projections of an object from dierent angles and then reconstruct the 2-D cross-sectional image through more or less complex computations. The basic mathematical principles behind the reconstruction of CT images was proposed by the Austrian mathematician Radon. In 1917, he proposed the famous Radon transform which transforms a function into its line integrals over particular lines, and also provided the inverse Radon transform which can be used to reconstruct images from a number of projections [START_REF] Radon | Uber die bestimmung von funktionen durch ihre integralwerte langs gewissez mannigfaltigheiten[END_REF][START_REF] Radon | On the determination of functions from their integral values along certain manifolds[END_REF]. In 1963, a South African American physicist

Cormack introduced a method to nd the real function in a nite region of a plane given its line integrals along all straight lines intersecting the region, he also pointed out that this method can be used to determine 2 dimensional X-ray attenuation coecients [Cormack, 1963]. However, he never considered of putting his theory into practice.

During the 1960s and early 1970s, there were a lot of researchers who dedicated on the development of CT technology [Kalender, 2006], including Oldendorf [Oldendorf, 1961],

Kuhl and Edwards [START_REF] Kuhl | Image separation radioisotope scanning[END_REF], Bracewell and Riddle [START_REF] Bracewell | Inversion of fanbeam scans in radio astronomy[END_REF], Gordon [START_REF] Gordon | Algebraic reconstruction techniques (art) for three-dimensional electron microscopy and x-ray photography[END_REF], Bates and Peters [START_REF] Bates | Towards improvements in tomography[END_REF]]. The most important breakthrough was made by an English electrical engineer Hounseld. In 1971, he invented the rst X-ray CT scanner, the EMI Mark I, and installed it at the Atkinson Morleys Hospital in Wimbledo (see Figure 2.9 (a)). The rst patient examination using this CT scanner took place in the same year and the obtained images revealed the presence of cyst formation (Figure 2.9 (b)) [Hounseld, 1973, Ambrose, 1973].

The great success of the rst head CT scanner instantly attracted numerous neurologists, radiologists, physicists, engineers and data processing specialists, they started to work on the algorithms of CT reconstruction and the interpretation of CT images. Commercial manufacturers like General Electric (GE), Siemens and Philips soon joint in the race to seize more market for this fresh but prospective technology. The number of CT scanners increased from 60 in 1974 to more than 10000 in 1980 [Kalender, 2006]. To reward the great achievements of Hounseld, he was awarded the Nobel Prize for Physiology or Medicine in 1979 along with Cormack.

Afterwards, the advancement of X-ray CT has greatly beneted from developed recon- The introduction of spiral (helical) CT in 1990 [START_REF] Kalender | Spiral volumetric ct with single-breath-hold technique, continuous transport, and continuous scanner rotation[END_REF] is another important step in the development of CT scanners. Spiral scan mode made it possible to continuously scan the extended volumes without interscan delay. The advantages of speed and continuity oered improved CT applications, such as the imaging of the lung within one breathhold. In 1998, multi-slice CT was proposed [START_REF] Klingenbeck-Regn | Subsecond multi-slice computed tomography: basics and applications[END_REF], Hu, 1999], even larger volume coverage in shorter scan times and improved longitudinal resolution became feasible [START_REF] Ulzheimer | Multislice ct: current technology and future developments[END_REF]. Later on, the 8-slice [START_REF] Kachelrieÿ | Ecgcorrelated image reconstruction from subsecond multi-slice spiral ct scans of the heart[END_REF], 16-slice [START_REF] Flohr | New technical developments in multislice ctpart 1: Approaching isotropic resolution with sub-millimeter 16-slice scanning[END_REF], Flohr et al., 2002a], 64-slice (2004, GE, Philips and Toshiba), 256-slice and 320-slice (2007, Philips and Toshiba) spiral CT scanner were implemented successively. It is striking to see how much the CT technology has evolved since its rst introduction several decades ago.

X-ray CT was rstly invented for medical application, and it is now an important medical imaging tool for preventive medicine and the display of lesions. For example, a cardiac CT scan for calcium scoring gives pictures of the coronary arteries to determine plaques, or contrast (usually iodine) enhanced heart scan to assess the extent of occlusion in the coronary arteries. CT is also used in the applications outside medical area. For example, defects detection in lumber [START_REF] Funt | Detection of internal log defects by automatic interpretation of computer tomography images[END_REF], exploration of coal [START_REF] Mathews | A review of the application of x-ray computed tomography to the study of coal[END_REF], inspection of ber reinforced polymer [START_REF] Schilling | X-ray computed microtomography of internal damage in ber reinforced polymer matrix composites[END_REF], Kastner et al., 2010] and geology engineering [START_REF] Mees | Applications of x-ray computed tomography in the geosciences[END_REF]. For convenience, the term CT in this thesis is used solely to refer to X-ray computed tomography.

X-ray source

X-ray can be generated through the Bremsstrahlung process. We can accelerate the electrons in a vacuum tube with a high voltage U acc , and then stop the high velocity electron with a metal target (usually tungsten, alloy of rhenium and tungsten, molybdenum or copper). The deceleration of electrons will lead to the emission of Bremsstrahlung, where the lost kinetic energy is converted into photons with continuous spectrum. The generated photons can be in any directions and their maximum energy will be E max = eU acc , the kinetic energy of electrons. For example, if the tube voltage is 100 kV, the generated X-rays will have energy up to 100 keV. In addition, the low-energy photons have weak penetration, they are likely to be fully absorbed by large object (like human) and can not arrive the detector, thus they contribute to patient dose but not to the image production. Therefore, the low-energy photons are usually removed by a lter in front of the X-ray tube before irradiating any object [Walsh, 2014].

On the other hand, the incident electrons can also interact with the orbital electrons.

When an orbital electron is knocked out of the inner shell of the target metal atom, the electron from higher energy levels will ll the inner shell vacancy, resulting in emission of X-ray. Usually the transitions are from upper shells into K shell (called K lines), into L shell (called L lines). The wavelength of X-rays emitted through this process depends on the type of metal material and the energy gap between the respective electron orbital shells. Thus the X-rays have several specic discrete energies, called characteristic X-ray.

Figure 2.10 shows the X-rays generated from a tungsten anode tube with voltage of 150 kV, the Bremsstrahlung emission and characteristic X-rays can be observed in this spectrum. [Cattin, 2016] The intensity of emitted X-rays are determined by photon ux and the energy of photons. The photon ux is dened as the number of photons per second per unit area and it is proportional to the tube current. The photon counting process has Poisson distribution, that is the variance of photon counting equals the expected counting. The signal-to-noise ratio dened by the ratio of mean to standard deviation of the counting is proportional to the square root of the expected counting. Thus increasing the photon ux will gain more signal quality.

Energy integrating detector

Conventional CT devices are equipped with energy integrating detectors which provide integrated information of all the transmitted beam. Two types of energy integrating de-tectors prevail: scintillator-based detector and semiconductor-based detector.

Scintillator-based detector consists of three essential layers: scintillator layer, photodiode layer and a substrate to provide the mechanical and electrical infrastructure [START_REF] Shefer | State of the art of ct detectors and sources: a literature review[END_REF]. The physical principles are illustrated in 2.11. When the X-rays interact with the scintillator, photons are absorbed predominantly via the photoelectric eect, and a number of electron-hole pairs are created. After some time the recombination of electrons with holes will result in the emission of photons within visible range [Nikl, 2006].

The

CT image reconstruction

Image reconstruction is an important step to retrieve the cross sectional attenuation coecients from the measured X-ray data, which is the transmitted beam intensity along a line between X-ray source and detector. For monochromatic X-ray source of energy E, the measurement can be written as a function of the attenuation coecients: (2.14) where I(E) is the measured intensity and I 0 (E) is the intensity coming from the X-ray source, µ E (x, y) is the linear attenuation coecient of the object at location (x, y). Figure 2.13 represents the case where a parallel beam is used, in which θ and t are polar coordinates of the path that the X-ray beam crosses through the object, like in Figure 2.13 is. Thus the line integration of the µ E (x, y) is:

I(E) = I 0 (E) exp(- ∫ L µ E (x, y) ds),
P E,θ (t) = ∫ (θ,t)line µ E (x, y) ds = -ln I(E) I 0 (E) .
(2.15)

In practical cases, the X-rays are polychromatic. However, the Equation (2.15) is also applied and this yields artifacts in the CT reconstruction due to beam hardening since we usually ignore the inuence of dierent photon energies. For the following part in this section, we will omit the energy dependent variable E, since it does not contribute to the method of reconstruction.

We generallize µ(x, y) with f (x, y), the path that the X-ray crosses in the object can be given by:

x cos θ + y sin θ = t.

(2.16)

Thus Equation (2.15) can be rewritten as: To reconstruct µ(x, y) from P θ (t), there are two main categories of methods, analytical reconstruction and iterative reconstruction. The principle are reviewed in the book of [START_REF] Kak | Principles of computerized tomographic imaging[END_REF]. The most common algorithm in the analytical reconstruction category is the ltered back-projection (FBP) method.

P θ (t) = ∫ +∞ -∞ ∫ +∞ -∞ f (x, y)δ(x cos θ + y sin θ -t)

Filtered back-projection reconstruction method

We rotate the original coordinate system by angle θ, can change (x, y) into (t, s), the transform can be expressed as:

[ t s ] = [ cos θ sin θ -sin θ cos θ ] [ x y ] (2.18)
In the (t, s) coordinate system, we do a Fourier transform to P θ (t):

S θ (w) = ∫ +∞ -∞ P θ (t)e -j2πwt dt = ∫ +∞ -∞ [∫ +∞ -∞ f (t, s)ds ] e -j2πwt dt. (2.19)
According to Equation (2.18), the above formula can be written as:

S θ (w) = ∫ +∞ -∞ ∫ +∞ -∞
f (x, y)e -j2πw(x cos θ+y sin θ) dx dy.

(2.20)

Make u = w cos θ, v = w sin θ, then the right-hand side of Equation ( 2.20) is a twodimensional Fourier transform at spatial frequency of (u, v):

S θ (w) = F (w cos θ, w sin θ) = F (u, v).
(2.21)

The above equation indicates that the Fourier transform of the original object on radial lines can be determined by the Fourier transform of the projections at dierent angles. This is called the central or Fourier slice theorem ( or projection-slice theorem). Fourier slice theorem is the basis of FBP algorithm. From Equation (2.20) and (2.21), we know that: ux+vy) du dv. (2.24) equation ( 2.23) can be written as:

f (x, y) = ∫ +∞ -∞ ∫ +∞ -∞ F (u, v)e j2π(
f (x, y) = ∫ π 0 [∫ +∞ -∞ F (w, θ)|w|e j2π(x cos θ+y sin θ) dw ] dθ.
(2.25)

We substitute the two-dimensional Fourier transform F (w, θ) of object with the Fourier transform of the projection data at angle θ, and the above equation can be expressed as:

f (x, y) = ∫ π 0 [∫ +∞ -∞ S θ (w)|w|e j2π(x cos θ+y sin θ) dw ] dθ = ∫ π 0 Q θ (x cos θ + y sin θ) dθ, (2.26)
where

Q θ (x cos θ + y sin θ) = Q θ (t) = ∫ +∞ -∞
S θ (w)|w|e j2πt dw.

(2.27)

Equation (2.27) indicates that the original projection data P θ (t) is rstly transformed into frequency domain S θ (w), then S θ (w) is multiplied by lter |w|, and the result is transformed into Q θ (t) by the inverse Fourier transform. Therefore Q θ (t) is considered as ltered projections". Equation (2.27) indicates that the ltered projections Q θ (t) are back-projected". Through integration of θ, the ltered projections of all directions are smeared back for the determination f (x, y).

The ltering step can also be done in space domain using convolution:

Q θ (t) = ∫ +∞ -∞ P θ (t ′ )h(t -t ′ ) dt ′ , (2.28)
where kernel h(t) is the inverse Fourier transform of |w|. In fact, the inverse transform of |w| does not exist since it is not square integrable. Since the projection data is sampled, the lter only needs to match |w| at the sampling points. Also, in practice the energy contained in high frequency is negligible, the lter can be made bandlimited. Suppose the W is the maximum frequency, and according to the sampling theorem, the sampling interval is τ = 1/2W . Bandlimited lter can be expressed as:

H(w) = |w|b w (w), (2.29) b w (w) = { 1 |w| < W 0 otherwise.
(2.30)

The corresponding lter kernel of H(w) is:

h(nτ ) =      1 4τ 2 , n = 0 0, n even -1 n 2 π 2 τ 2 , n odd.
(2.31) Function H(w) and function h(nτ ) are shown in Figure 2.15. In practice, a modied lter will be used to weaken the higher frequencies and obtain better image quality. This can be achieved by a simple multiplicative window such as a Hamming window.

-1/2 0 1/2 Frequency w 0 1/2 H(w) (a) -5 -3 - 0 3 5 n 0 1/4 2 h(n ) -1/(25 2 2 ) -1/(9 2 2 ) -1/( 2 2 ) (b)
In practice, fan beam CT is more common than parallel CT. The X-ray source radiates fan-shaped beam and a number of detectors are put opposite to detect the transmitted photons. They rotate together to generate a set of projections. Two types of fan beam projections exist, like discussed in [START_REF] Kak | Principles of computerized tomographic imaging[END_REF]: the equiangular interval projections and the equispaced interval projections. The former applies for conditions where the detectors are put to sample the projections by equiangular intervals; and the latter is for the condition where the detectors are positioned with equal space along a straight line (see Figure 2.16). In our research, the spectral CT imaging system accords with the latter one, therefore we will discuss the FBP reconstruction algorithm in this situation. We assume that the projection data acquired from Figure 2.16 is R β (s), and its relationship with the parallel projections P θ (t) can be derived from Figure 2.17:

X-ray source

Equally spaced detectors

t = s cos γ = sD √ s 2 + D 2 (2.32) θ = β + γ = β + tan -1 s D .
(2.33)

Similarly with Equation 2.26, the reconstruction algorithm can be derived for fan beam CT imaging as:

R ′ β (s) = R β (s) D √ D 2 + s 2 (2.34) g(s) = 1 2 h(s) (2.35) Q β (s) = R ′ β (s) * g(s) (2.36) f (r, ϕ) = ∫ 2π 0 1 U 2 Q β (s ′ ) dβ, (2.37)
where r and ϕ are polar coordinates of a pixel on the object, s ′ is the value of s for the ray that passes through the pixel (r, ϕ) under consideration, and U is the ratio between distance D and the distance from the source to the projection of point (r, ϕ) onto the central ray:

U = D + r sin(β -ϕ) D (2.38) s ′ = D r cos(β -ϕ) D + r sin(β -ϕ)
.

(2.39) γ β θ X-ray source

Iterative reconstruction method

Another category of image reconstruction for X-ray CT is iterative reconstruction method.

This method begins with an initial image assumption, and compares it to measured values while making adjustments until the two are in agreement. by the following dierent steps:

1 The initial image f (0) is arbitrary, but usually we make it 0, or an image reconstructed by another method such as FBP. The stopping criterion can either be the maximum number of iterations, or the small enough update of current image estimate , or the predened quality criterion [START_REF] Beister | Iterative reconstruction methods in x-ray ct[END_REF]. A set of iterative reconstruction methods has been proposed, such as the algebraic reconstruction technique (ART) [START_REF] Gordon | Algebraic reconstruction techniques (art) for three-dimensional electron microscopy and x-ray photography[END_REF],

⃝
simultaneous iterative reconstruction technique(SIRT) [Gilbert, 1972], simultaneous algebraic reconstruction technique (SART) [START_REF] Andersen | Simultaneous algebraic reconstruction technique (sart): a superior implementation of the art algorithm[END_REF], maximum likelihood expectation-maximization (ML-EM) algorithm [START_REF] Shepp | Maximum likelihood reconstruction for emission tomography[END_REF], and so on. Here we introduce in detail the most basic and popular ART method.

In ART reconstruction, the projection process of CT imaging system can be represented by:

b i = ∑ j a i,j f j , (2.40)
where b i is the integration data measured with the i th ray (total number of rays M =No.rotation angles × No.detectors); a i,j is the weighting factor that represents the contribution of the j th voxel to ray i; f j is the attenuation coecient of voxel j in our case. Therefore, we have:

Af = b, (2.41) where A = a i,j is the system matrix, f = [f 1 , f 2 ..., f N ] is the image vector where N is the total number of voxels in the reconstruction image, and b = [b 1 , b 2 ..., b M ] is the measurement vector. To solve f , the Kaczmarz method [Kaczmarz, 1993] is applied and the following iterative formula is obtained:

f (k) = f (k-1) + λ k-1 b i -f (k-1) a T i a i a T i a i , (2.42)
where k is the number of iterations and λ k-1 is the relaxation parameter, i = k mod M + 1, a i is a vector made of the i th row of A, .

Compared to analytical FBP method, iterative reconstruction has advantage to reduce noise and improve image quality by integrating various physical models and incorporating prior information. It outperforms FBP method especially with sparse and incomplete data.

The most serious disadvantage of iterative reconstruction is its huge computational costs, since each iteration involves forward projections. Computing power may not be sucient for complicated modeling, such as in the case of ne volumetric image reconstruction, dualenergy or spectral CT data processing, extreme large matrix sizes will be used and great demands on computational speed will be requested. For this reason iterative reconstruction techniques are not yet suitable for scanning patients in acute situations [Opie, 2013].

2.3 X-ray spectral CT

Introduction

We have discussed in the above the principles of conventional CT imaging. The X-ray sources in use nowadays consists in a spectrum of energies. However, the conventional energy-integrating detector sums up all the energy information and thus gives lower energy photons which carries contrast information less weight than higher energy photons (Figure 2.18 (a)). The reconstructed image corresponds to the linear attenuation coecients at a single average energy. Spectral information is not fully utilized in this case and contrast between dierent materials is lost. Spectral CT appears to make advantage of multi-energy information passed through the object to improve CT imaging.

The concept of energy-resolving CT has been discussed since the invention of CT. In the landmark CT paper [Hounseld, 1973] of Hounseld, he proposed to take CT scans with the same slice of material at two dierent voltages to determine atomic number of materials. In 1976, Alvarez and Macovski proposed to extract energy-dependent information by representing the total attenuation with functions corresponding to photoelectric and Compton eect, and the energy-independent parts of the two eects were then reconstructed respectively [START_REF] Alvarez | Energy-selective reconstructions in x-ray computerised tomography[END_REF]. In 1977, Riederer et al. presented their work of selective iodine imaging by CT scans using three heavily ltered X-ray beams [START_REF] Riederer | Selective iodine imaging using k-edge energies in computerized x-ray tomography[END_REF]. Material selective imaging using dual-energy scan also extended to X-ray radiographic imaging [START_REF] Brody | Dual-energy projection radiography: initial clinical experience[END_REF], Lehmann et al., 1981].

Following the above concept, a variety of dual-energy CT (DECT) systems were built for energy-resolving imaging. There are currently three main techniques for achieving DECT. The rst one is called rapid kilovoltage switching technique, where the voltage of X-ray tube is switched rapidly to produce two X-ray beams with dierent energies for scanning at every other projection during the whole scan, and its principle is shown in and photodiodes stacked on top of each other, so that the lower-energy photons are absorbed by the top scintillator array, while the higher-energy photon penetrates the top layer and are absorbed by the bottom scintillator array [Nasirudin, 2015], see Figure 2.18 (c). DECT has limited ability in material decomposition by producing two datasets corresponding to two dierent X-ray spectra of the same anatomic region [START_REF] Danad | New applications of cardiac computed tomography: dual-energy, spectral, and molecular ct imaging[END_REF].

The potential of using DECT to detect the distribution of materials in vivo attracted much investigation [START_REF] Fischer | Quantication of liver iron content with ctadded value of dual-energy[END_REF], Hazirolan et al., 2008, Luo et al., 2015, Panta et al., 2015].

Take an essential element for human body iron as an example: The main advantage of PCD over energy integrating detectors is the eliminated electronic noise. By setting a low-energy threshold that is above the highest amplitude of the noise, electronic noise can be excluded from being counted. In this way, low-dose and low-energy benets with higher image quality. Furthermore, instead of being less weighted, low-energy photons can be counted in order to bring valuable information, e.g. more contrast to the image. However, PCDs are not awless, the measurement of PCDs may deviate from the truth for several reasons [START_REF] Taguchi | Vision 20/20: Single photon counting x-ray detectors in medical imaging[END_REF]]:

• Pulse pileup. When several photon interactions occur at the same time, multiple pulses may be piled up and be counted as one event, resulting in wrong recorded energy and counts loss.

• Charge sharing. When the generated electron charge cloud reaches the anode near the pixel boundary, they maybe divided and counted by the neighboring pixels to a smaller energy level.

• K-escape X-rays. X-ray uorescence or Auger electrons may be generated during the photoelectric eect, as we have discussed in 2.1.1. The uorescent X-ray photons can either be absorbed by the PCD pixel with the primary interaction again, be detected by an adjacent pixel, or leave the PCD completely.

• Compton scattering. The photon energy is not fully transmitted to electrons, resulting in energy loss.

• Charge trapping. The electron or hole is trapped by a trapping center and is thermally re-emitted with a delay. This will result in lower-energy counting.

• Polarization and long-term reliability. Reference group 5 [Fredenberg et al., 2010a[START_REF] Fredenberg | Energy resolution of a photon-counting silicon strip detector[END_REF] 6

Reference group 6 [START_REF] Ballabriga | Medipix3: A 64 k pixel detector readout chip working in single photon counting mode with improved spectrometric performance[END_REF], Ballabriga et al., 2006, Gimenez et al., 2011, Zainon, 2012] 7 Reference group 7 [START_REF] Alessio | [END_REF]MacDonald, 2013, Hamamatsu, 2014] Page 32

Cette

Quantitative material decomposition

The X-ray attenuation for dierent materials varies from each other. With the spectral information of transmitted photons, it is possible to discriminate dierent components of object. As we have said before, in 1976, Alvarez and Macovski proposed to extract energydependent information by represent the total attenuation with functions correspond to photoelectric and Compton eect, and the energy-independent parts of the two eects were then reconstructed [START_REF] Alvarez | Energy-selective reconstructions in x-ray computerised tomography[END_REF]. The reconstructed photoelectric image provides composition information of object while the reconstructed Compton image provides mass density information of object. In 2007, Roessl and Proksa developed this idea by adding the K-edge terms on the photoelectric and Compton eect decomposition formula and gadolinium was discriminated from normal tissues in the thorax phantom in their simulation study [START_REF] Roessl | K-edge imaging in x-ray computed tomography using multi-bin photon counting detectors[END_REF]. Their results were further validated by CT system experiments where the contrast agents (gadolinium and iodine) were well quantied [START_REF] Schlomka | Experimental feasibility of multi-energy photon-counting k-edge imaging in pre-clinical computed tomography[END_REF]. This technique was named K-edge imaging and it enables multi-agent imaging [START_REF] Roessl | Preclinical spectral computed tomography of gold nano-particles[END_REF], Butler et al., 2011, Boone et al., 2014, Pan et al., 2010, Si-Mohamed et al., 2017, Lee et al., 2016b, Lee et al., 2016a].

There is less research about the quantication ability of materials that lack K-edge properties within the detection energy range, which is more challenging. water-like, calcium-like and lipid-like tissues were well separated according to histology comparison, but iron was not distinguished from calcium [Zainon, 2012].

Detailed material decomposition methods will be introduced in Chapter 3.

Energy weighting

There are mainly three weighting cases in CT imaging: integrating, counting and energy weighting. The rst case uses weighting factors that are proportional to the energy of photon; the second case sums up the number of interacting photons with constant weighting factor; and in the third case, the energy-bin data are weighted by energy dependent factor to obtain better contrast-to-noise ratio (CNR) and signal-to-noise Ratio (SNR) [START_REF] Carramate | Energy weighting technique in quantum computed tomography using a mpgd[END_REF]. By giving lower-energy data more weight, the dose eciency is eectively improved and enables a higher CNR at a constant patient dose, or a lower dose at a constant CNR [START_REF] Berglund | Energy weighting improves dose eciency in clinical practice: implementation on a spectral photon-counting mammography system[END_REF]. Energy weighting can be performed both in projection domain or image domain data. In the research of Giersch el al., by using weighting factor E -3 , the SNR can be enhanced by up to from 1.3 to 1.9 or reduce the dose by factor of 2.5

without losing image quality [START_REF] Giersch | The inuence of energy weighting on x-ray imaging quality[END_REF]]. An image-based" weighting method by Schmidt shows that the CNR is improved by factors of 1.0 to 1.3 compared to photoncounting weighting with negligible beam-hardening artifacts [Schmidt, 2009].

Virtual monochromatic imaging

We have discussed in 2.1 the principles of X-ray interactions with matter. There are limited number of interaction eects and their energy dependence have all been represented.

Therefore through the energy-resolving imaging, virtual monochromatic images can be calculated in the projection domain or in the image domain. Such monochromatic images have the advantage to improve spectral CNR, since the low-energy attenuation carries more contrast. For example, in intra-hepatic and extra-hepatic portal veins imaging, 51 keV monochromatic images provide the best CNR that is 100% higher than the polychromatic images [START_REF] Zhao | Improving image quality in portal venography with spectral ct imaging[END_REF]. Hard X-rays attenuate less than soft X-rays, thus monochromatic images of higher energy have the potential to reduce beam-hardening [START_REF] Neuhaus | Improvement of image quality in unenhanced dual-layer ct of the head using virtual monoenergetic images compared with polyenergetic single-energy ct[END_REF], Yu et al., 2012] and metal artifacts [START_REF] Wang | Metal artifacts reduction using monochromatic images from spectral ct: evaluation of pedicle screws in patients with scoliosis[END_REF]. Another application of monochromatic imaging is to characterize materials. For example, when comparing a set of monochromatic images at lower energies, the iodine and calcium can be separated since the attenuation of iodine changes greatly around its K-edge (33.17 kev) while that of calcium is quite at. Therefore the CNR of contrast-enhanced imaging can be improved by choosing an energy slightly above the absorption edge of the contrast agent.

Summary

This chapter presents the general principles of spectral CT, including the following sections:

1. Interactions of X-rays with matter. X-ray is a form of electromagnetic radiation that has extremely short wavelength and high frequency. It was discovered in 1985 by Röntgen. The research on new imaging modality using X-rays started since its discovery. Those research are based on the attenuation property of X-rays when passes through the object. There are three major X-ray interactions with matters that account for the total attenuation: photoelectric absorption, Compton scattering and pair production. Their mechanisms are introduced in detail. The rst is eect-based decomposition, also called K-edge imaging, where the attenuation is considered as a linear combination of the photoelectric eect, Compton eect and extended K-edge components if there are materials with their K-edge within the CT energy range. The second is material-based method, which describes the linear attenuation coecient by the contributions of a basis of materials.

3.1.1 Eect-based decomposition / K-edge imaging K-edge is the binding energy of the K shell electron of an atom. Generally, the attenuation coecient decreases with photon energy, but there is a sudden increase at the photon energy just above the binding energy since the probability of the photoelectric eect occurring is much higher than the energy below this binding energy. There exists also L-edge which corresponds to the binding energy of L shell electrons. Figure 3.1 illustrates the energy dependence of X-ray attenuation coecient for iron (a) and iodine (b), where can be clearly seen the K-edges of the two materials as well as the L-edges of iodine. A K-shell interaction is approximately four to ve times more probable than an L-shell interaction if both interactions are energetically allowed, thus we consider only the K-edge here. It should be notice that the K-edge of iron appears at a quite small energy (7.11 keV) that is usually outside the energy range of PCD, therefore we dene such lighter materials that lack K-edge properties within the detection energy range as non-K-edge material", and those having K-edge absorption within the energy range as K-edge material". Generally, K-edge materials are much easier to be quantied because of their unique attenuation coecient. K-edge of iron appears at 7.11 keV; iodine K-edge at 33.17 keV. L-edges of iodine : L-I 5.19 keV, L-II 4.85 keV and L-III 4.56 keV. Data taken from [START_REF] Hubbell | Tables of x-ray mass attenuation coecients and mass energy-absorption coecients (version 1.4)[END_REF].

The detector absorption eciency is denoted as D(E) and we assume that the photon counting detector allows energy-resolving measurements in N energy bins with ideal energy resolution. According to Beer-Lambert law, the expected number of photons λ i , in energy bin B i (i=1, 2, ..., N ) can be expressed as

λ i = E f (i) ∑ E=Es(i) D(E)N 0 (E) exp[- ∫ µ( -→ x , E) ds], (3.1)
where E s (i) and E f (i) denote the start and nal energies of bin B i , N 0 (E) is the number of photons in the initial spectrum at energy E, µ( -→ x , E) represents linear attenuation coecient of the object at location -→ x and energy E. It is noticed that there is a location-dependent integral and an energy-dependent sum in this formula, which correspond respectively to the two variables -→ x and E on which µ( -→ x , E) depends. Therefore, if we transform µ( -→ x , E) into the sum of the products of energy-dependent term f α (E) and location-dependent term a α ( -→ x ), it will be easier to calculate the integral and the sum separately.We then write

µ( -→ x , E) = M ∑ α=1 a α ( -→ x )f α (E), (3.2)
where M denotes the number of terms in use. Generally, a α ( -→ x ) describes the physicochemical or quantity property (depending on the concrete formulation of Equation (3.2)) of the unknown object and f α (E) is a known function that varies with energy. In the work of [START_REF] Roessl | K-edge imaging in x-ray computed tomography using multi-bin photon counting detectors[END_REF], µ( -→ x , E) is taken as the combination of photoelectric absorption, Compton scattering and extra K-edge components if there is any K-edge material in the object, making Eq. ( 3.

2) into µ( -→ x , E) = a ph ( -→ x )f ph (E) + a Co ( -→ x )f KN (E) + M -2 ∑ α=1 ρ α ( -→ x )µ mα (E), (3.3)
where a ph ( -→ x ) and a Co ( -→ x ) represent location dependency of photoelectric absorption and Compton scattering, respectively. f ph (E) denotes the cross section of photoelectric absorption that can be approximated by E -3 , see Equation (2.6). f KN (E) is the cross section of Compton scattering, namely Klein-Nishina function, see Equation (2.10). ρ α ( -→ x ) denotes the density of K-edge material α at point -→ x and µ mα (E) the mass attenuation coecient of K-edge material α at energy E. (M -2) represents the number of K-edge materials inside the object.

Substituting Equation (3.2) into Equation (3.1), we have

λ i (A α ) = E f (i) ∑ E=Es(i) D(E)N 0 (E) exp[- ∫ M ∑ α=1 a α ( -→ x )f α (E) ds] = E f (i) ∑ E=Es(i) D(E)N 0 (E) exp[- M ∑ α=1 f α (E) ∫ a α ( -→ x ) ds] = E f (i) ∑ E=Es(i) D(E)N 0 (E) exp[- M ∑ α=1 f α (E)A α (sx, sy)], (3.4) with A α (sx, sy) = ∫ a α ( -→ x ) ds, (3.5)
where A α (sx, sy) is the line integral of location-dependent parameter a α ( -→ x ) along the measured projection path, (sx, sy) represents the location of the measured projection in the sinogram, with sx denoting detector pixels and sy denoting measurement angles. To estimate A α (sx, sy), maximum-likelihood (ML) method is used in [Roessl andProksa, 2007, Schlomka et al., 2008]. Firstly suppose the detected number of photons M i form a set of independent Poisson random variables,then the likelihood function as the possibility of measurement result (M 1 = m 1 , M 2 = m 2 , ... M N = m N , ) can be calculated: Aα) .

P (m 1 , ..., m N |λ 1 (A α ), ..., λ N (A α )) = N ∏ i=1 λ i (A α ) m i m i ! e -λ i (
(3.6)

For computational convenience in seeking the extremum, Equation (3.6) is simplied by applying a negative logarithm operator on both sides and dropping the term that is unrelated to A α . Then calculating the negative logarithm likelihood function, we will obtain:

L(A α ) = -ln [P (m 1 , ..., m N |λ 1 (A α ), ..., λ N (A α )] = N ∑ i=1 [λ i (A α ) + ln m i ! -m i ln λ i (A α )] ∼ = N ∑ i=1 [λ i (A α ) -m i ln λ i (A α )] .

Material-based decomposition

There is actually another commonly used standard to categorize the material decomposition methods for spectral CT. According to the type of data on which the decomposition step operates, the methods can be divided into projection domain method and image domain method. The former, also called pre-reconstruction decomposition method, decomposes projection data into several desired components (i.e. materials), and then applies independently reconstruction to each component to obtain their respective distribution in space. The latter, also known as post-reconstruction method, operates on conventionally reconstructed CT images corresponding to dierent energy bins to obtain the concentration of target materials using linear matrix inversion, see gure 3.2. Both methods have their advantages and disadvantages, depending on dierent situations. The K-edge imaging presented in 3.1.1 is actually a projection domain method since it processes directly the projection data. We know from K-edge imaging that only the K-edge materials with relatively high atomic number can be discriminated using the projection domain method. However, in practice, non-K-edge materials are more commonly encountered. Lighter atoms, such as iron, calcium and potassium, usually have much lower K-edge energy that is beyond the detection energy range. In this situation, K-edge imaging method can not distinguish them, nor give any quantitative information since their attenuations are all included and mixed in photoelectric absorption and Compton scattering. Therefore, another version of formula Equation (3.2) should be considered to describe linear attenuation coecient in terms of a basis of M materials. We will develop the material-based method in two situations: projection domain decomposition and image domain decomposition.

Projection domain material-based method

We have introduced that the attenuation coecient can be considered as the linear combination of the basis material mass attenuation coecients weighted by their densities. So

we have µ( -→ x , E) = ρ 1 ( -→ x )µ m1 (E) + ρ 2 ( -→ x )µ m2 (E) + ... + ρ M ( -→ x )µ mM (E) = M ∑ α=1 ρ α ( -→ x )µ mα (E), (3.8)
where ρ α ( -→ x ) denotes the density of material α at point -→ x . Material α is not exclusive to K-edge materials, but also applicable for lighter materials. µ mα (E) designates the mass attenuation coecient of material α at energy E, which can be found in NIST database [START_REF] Hubbell | Tables of x-ray mass attenuation coecients and mass energy-absorption coecients (version 1.4)[END_REF]]. M should be smaller than or equal to N , otherwise the problem will be under-determined and an innite number of solutions exist if there are no other equality constraints. The types of materials chosen for decomposition depend on the application and prior knowledge of the scanned objects.

Similarly with Equation (3.4), by substituting Equation (3.8) into Eq. (3.1), we have

λ i = E f (i) ∑ E=Es(i) D(E)N 0 (E) exp[- ∫ M ∑ α=1 ρ α ( -→ x )µ mα (E) ds] = E f (i) ∑ E=Es(i) D(E)N 0 (E) exp[- M ∑ α=1 µ mα (E) ∫ ρ α ( -→ x ) ds] = E f (i) ∑ E=Es(i) D(E)N 0 (E) exp[- M ∑ α=1 µ mα (E)P α (sx, sy)],
(3.9)

with 3.10) where P α (sx, sy) is the line integral of object density ρ α ( -→ x ) along the measured projection path. Normally the objective function used to estimate P α (sx, sy) is the least squares (LS) method: (3.11) where m i is the measured number of photons in energy bin B i . After P α is obtained, we can apply conventional reconstructions to get the density distribution of target materials.

P α (sx, sy) = ∫ ρ α ( -→ x ) ds, (
P α = arg min Pα N ∑ i=1 [λ i (P α ) -m i ] 2 ,
We denote this method as proj-LS method.

Image domain material-based method

Unlike projection domain method, decomposition in image domain works on conventionally reconstructed images µ(i, -→ x ) for each energy bin B i , which represents energy-resolved linear attenuation coecient. Similarly to Equation (3.8), µ(i, -→ x ) is considered to be the linear combination of mass attenuation coecients µ mα (i) within one energy bin B i weighted by the concentration of material ρ α ( -→ x ):

µ( -→ x , i) = ρ 1 ( -→ x )µ m1 (i) + ρ 2 ( -→ x )µ m2 (i) + ... + ρ M ( -→ x )µ mM (i) = M ∑ α=1 ρ α ( -→ x )µ mα (i).
(3.12)

The main dierence between Equation (3.12) and Equation (3.8) is that µ mα (i) denotes the mass attenuation coecient of material α within energy bin B i , instead of a punctual xed energy value. It is not trivial to determine the exact µ mα (i) values since precisely it correlates not only with energy, but also with unknown object properties. A simplied formula to determine the eective mass attenuation coecient of one energy bin from punctual energy values is dened by [START_REF] Le | Least squares parameter estimation methods for material decomposition with energy discriminating detectors[END_REF]:

µ mα (i) = ∑ E f (i) E=Es(i) N 0 (E)µ mα (E) ∑ E f (i) E=Es(i) N 0 (E)
.

(3.13) By writing the location vector -→ x into two-dimensional coordinates (x, y), the density distribution of material α can be expressed as ρ α (x, y). And the least squares minimization to estimate ρ α (x, y) can be written as (3.14) where μ(x, y, i) represents the estimated linear attenuation coecient corresponding to energy bin B i from conventionally reconstructed images. This method is in image domain and is based on LS criterion, thus we call it ima-LS method.

ρ α (x, y) = arg min ρα N ∑ i=1 [μ(x, y, i) - M ∑ α=1 ρ α ( -→ x )µ mα (i)] 2 ,
Material-based decomposition has wider applications and also presents the potential to be optimized. Le and Molloi used a calibrated least squares tting technique in image domain to quantify calcium and iodine [START_REF] Le | Least squares parameter estimation methods for material decomposition with energy discriminating detectors[END_REF] for breast imaging. To assess tissue composition of atherosclerotic plaque, Alessioand and MacDonald used weighted least-squares method, which assumes that the materials are limited at a certain location to only those contained in one of the pre-selected classes. Therefore, this method rstly segments the image into pre-selected classes, then decomposes the pixel into corresponding materials contained in the class [START_REF] Alessio | Quantitative material characterization from multi-energy photon counting ct[END_REF]. Zeng et al. proposed a reconstruction method that uses penalized weighted least-squares scheme incorporating the structure tensor total variation regularization to gain better material decomposition quality [START_REF] Zeng | Penalized weighted least-squares approach for multienergy computed tomography image reconstruction via structure tensor total variation regularization[END_REF]. Ducros et al. uses a regularized weighted least squares Gauss-Newton algorithm to decompose the spectral X-ray projection images of a thorax phantom into soft tissue, bone and gadolinium bases [START_REF] Ducros | Regularization of nonlinear decomposition of spectral x-ray projection images[END_REF].

In the present work, material-based method is preferred according to our application. The commonly used objective function which is based on LS criterion (Equation (3.11) and (3.14)) directly minimizes the errors between expected and measured number of photons or the errors between expected and reconstructed attenuation coecients. However, the inverse problem of material decomposition is usually ill-posed and the X-ray spectral CT measurements suer from Poisson photon counting noises. The standard least squares objective function can lead to overtting to the noisy measurement data. We then propose to add a regularization term to penalize the sum of the square variations of the decomposed images, in order to reduce the eect of noise and enforce smoothness.

Projection domain

For projection domain decomposition, Equation (3.11) gives the objective function to solve the density integrals of dierent materials by minimizing the square errors between expected number of photons λ i and measured number of photons m i for all energy bins. According to Beer-Lambert law, the number of transmitted photons is negatively exponential to the attenuation coecient. Instead of minimizing directly the dierence of number of photons, we propose a least log-squares (LLS) criterion, as given in Equation (3.15). A natural logarithm operator is placed on λ i and m i , which makes them inversely related to the linear attenuation coecient. Therefore Equation ( 3.15) has a physical meaning of minimizing the errors on linear attenuation coecient. (3.15)

Material decomposition using the above LLS criterion is called proj-LLS method. Based on proj-LLS method, we further propose a patchwise regularized least log-squares (proj-PR-LLS) method. Firstly, the sinograms from acquisition are divided into patches with the size of np × np. The decomposition is performed patch by patch by minimizing the following objective function: 3.16) where P C α (sx, sy) represents P α (sx, sy) when (sx, sy) is inside patch C. r denotes the relaxation parameter. R(P C α ) is the regularization term used to reduce the eect of noise and enforce smoothness.

P C α (sx, sy) = arg min P C α (sx,sy) { ∑ (sx,sy)∈C N ∑ i=1 [ln(λ i (P C α )) -ln(m C i )] 2 + rR(P C α )}, (

R(p C

α ) in this work is considered as the sum of the L2 regularizations of gradient images of P C α :

R(P

C α ) = M ∑ α=1 || ▽ P C α || 2 2 = M ∑ α=1 ∑ (sx,sy)∈C (sx-1,sy)∈C (sx,sy-1)∈C {[P C α (sx, sy) -P C α (sx -1, sy)] 2 +[P C α (sx, sy) -P C α (sx, sy -1)] 2 }.
(3.17)

This regularization term enforces patchwise smoothness and can reduce the eect of noise.

Image domain

We have introduced in 3.1.2 the commonly used image domain material decomposition method, which can be summarized into two steps: rstly, reconstruct conventional (spatial) images of linear attenuation coecients from the acquired projection data corresponding to to dierent energy bins; then perform material decomposition on the reconstructed images.

Equation (3.14) gives a standard least squares objective function minimizing the square errors between expected and reconstructed linear attenuation coecients. To improve image quality, we propose to use the patchwise regularized least-squares (ima-PR-LS) method for the estimation of ρ α (x, y): (3.18) where ρ C α (x, y) is the density of material α when its location (x, y) is inside patch C. r denotes the relaxation parameter. R(ρ C α ) is the regularization term. Similarly to Equation 3.17, R(ρ C α ) is the sum of the L2 regularizations of gradient images of material density ρ C α :

ρ C α (x, y) = arg min ρ C α (x,y) { ∑ (x,y)∈C N ∑ i=1 [μ C (x, y, i) - M ∑ α=1 ρ C α (x, y)µ mα (i)] 2 + rR(ρ C α )},

R(ρ

C α ) = M ∑ α=1 || ▽ ρ C α || 2 2 = M ∑ α=1 ∑ (x,y)∈C (x-1,y)∈C (x,y-1)∈C {[ρ C α (x, y) -ρ C α (x -1, y)] 2 +[ρ C α (x, y) -ρ C α (x, y -1)] 2 }. (3.19)
Formula (3.18) minimizes the square errors between reconstructed attenuation coecient value and expected coecient value given by Equation (3.12), and at the same time constrains the total sum of squared variations within patches. In this chapter, we have described in detail the existing and proposed material decomposition methods for X-ray spectral CT, as illustrated in Figure 3.5.

According to the formulation of the decomposition, the existing material decomposition methods can be divided into two categories. The rst is eect-based decomposition, also called K-edge imaging, where the attenuation is considered as a linear combination of the photoelectric eect, Compton eect and extended K-edge components if there are materials with their K-edge within the CT energy range. The second is material-based method, which views the linear attenuation coecient as the contributions of a basis of materials. For the material-based method, we have introduced the least squares based decomposition in projection domain (proj-LS method) and in image domain (ima-LS method). The former method decomposes projection data into several desired components, and then applies independent reconstruction to each component to obtain the distribution of materials.

The latter operates on conventionally reconstructed CT images corresponding to dierent energy bins to obtain the concentration of target materials using linear matrix inversion.

Based on the existing methods, we have investigated the material-based approaches.

We have proposed a least log-squares criterion for projection domain method to minimize the errors on linear attenuation coecient: proj-LLS method. Furthermore, to reduce the eect of noise and enforce smoothness, we proposed patchwise regularized decomposition methods in both projection domain and image domain. The proposed methods minimize the square errors between the logarithmic expected and measured number of photons (proj-PR-LLS method), or between the reconstructed attenuation coecient value and expected coecient value (ima-PR-LS method), and at the same time constrains the total sum of squared variations within patches. The reason for having applied regularization to small patches instead of the whole image is to avoid huge computation cost. However, the

Introduction

The objective of this chapter is to introduce the simulation method of X-ray imaging and to evaluate the proposed projection domain material decomposition methods through the simulation study. Around these two points, rstly we will present the spectral CT simulation method based on Virtual X-ray imaging (VXI) software. Secondly, targeting on a medical application (atherosclerosis imaging, CT) and an industrial application (plastic sorting, computed radiography (CR)), we will simulate through VXI software the X-ray imaging process with two corresponding phantoms, apply the proposed methods and analyze their performance. For the former application, we will compare the decomposition performance of the proj-LS, proj-LLS and proj-PR-LLS method and for the latter application, we will use proj-PR-LLS method and a further developed proj-PR-LLS-clas method.

For both applications in this chapter, the spectral CT or CR imaging process are simulated with currently used acquisition parameters, e.g. number of energy bins and exposure factor. We will compare the performance of proj-PR-LLS and ima-PR-LS method, and investigate the inuence of these parameters on their performance in Chapter 5.

Simulation of spectral X-ray imaging 4.2.1 Basic principles

We simulate the spectral CT/CR system which has multiple energy bin resolving capability using Virtual X-ray imaging (VXI) software [START_REF] Duvauchelle | A computer code to simulate x-ray imaging techniques[END_REF]. VXI software was developed at INSA Lyon (France) to simulate the radiographic, radioscopic and tomographic imaging. The simulation is based on ray-tracing techniques together with the X-ray attenuation law.

N (E) = N 0 (E)∆Ω exp [ ∑ α -µ α (E)x α ] , (4.1)
where N 0 (E) is the emitted number of photons with energy E from the X-ray source per solid angle unit; ∆Ω is the corresponding solid angle of ray SK, which is proportional to the size of pixel K; µ α (E) represents the linear attenuation coecient of material α at energy E and x α the total path length through the material α. Photons scattered inside the object or the detector are not taken into account in the present simulation study.

VXI software allows the users to dene a set of parameters for their own applications.

X-ray source generation

In this software, multiple properties of X-ray source can be dened. Our simulation was carried out with point X-ray sources based on Birch & Marshall model [START_REF] Birch | Computation of bremsstrahlung x-ray spectra and comparison with spectra measured with a ge (li) detector[END_REF]. Figure 4.3 shows two spectra used in the present study: both were computed for tungsten target material without ltration, the target angles were 17 • and the tube voltage were respectively 100kVp and 120kVp. The beam shapes were chosen to be isotropic.

Other parameters like position and tube current were set respectively according to specic applications. Take water (H 2 O) as an example, we can dene this material either by assigning atomic percentages 66.67% and 33.33% or mass percentages 11.19% and 88.81% to hydrogen (H) and oxygen (O). Figure 4.4 illustrates the entire sample denition process.

Detector denition

Detector geometry can be determined by dening the number and size of pixels of the detector. It is also possible to dene the detector material (chemical composition, density, thickness) as well as the position and orientation of the detector plane. Targeted on the simulation of photon counting detector, we choose the spectro mode, and characterize the number of detector energy bins and the starting/nal energy of detecting energy range in order to allows counting the arriving photons into corresponding pre-dened energy bins to simulate the spectral X-ray imaging process.

Noise condition

During our simulation, photon noise was simulated for spectral CT/CR imaging. We will discuss in the following content the photon noise for polychromatic beam in photon counting mode, how we simulate it for spectral X-ray imaging and the inuence of the number of photons on the signal-to-noise ratio.

For monochromatic beam of energy E, the raw data coming from photon-counting process has Poisson distribution [Ochi, 1990[START_REF] Hasino | Photon, poisson noise[END_REF]. Denoting the expected number of photons as λ E and the measurement as m E , then the probability of detecting m E photons can be calculated by

m E ∼ Poisson(λ E ), p(m E ) = e -λ E λ m E E m E ! . (4.2)
Poisson distribution has the property that its variance σ 2 E equals the expected counting λ E , as given by σ 2 E = λ E . (4.3) This indicates that the photon noise depends on the expected number of photons λ E and that its standard deviation increases with the square root of λ E . When it comes to polychromatic beam, as in the spectral CT/CR imaging process, the detected number of photons m i for one energy bin B i is a sum of all the photons in that energy bin range (4.4) where E s (i) and E f (i) are the starting and nal energies of bin B i . Since m E (E = E s (i), ...E f (i)) are independent variables that obey the Poisson distribution, the expected value λ i and variance σ 2 i can be calculated

m i = E=E f (i) ∑ E=Es(i) m E ,
λ i = E=E f (i) ∑ E=Es(i) λ E , (4.5) σ 2 i = E=E f (i) ∑ E=Es(i) σ 2 E = E=E f (i) ∑ E=Es(i) λ E . (4.6)
It can be noticed that m i also has Poisson distribution: m i ∼ P oisson(λ i ). When λ i is greater than about 20, the Poisson distribution can be approximated by the Gaussian distribution with the mean and variance of λ i , and the approximation improves as λ i increases. In this case, we have

m i ∼ N (λ i , λ i ). (4.7)
Therefore, we can simulate the photon-counting process by adding a Gaussian noise to the expectation:

m i = λ i + σ i G = λ i + √ λ i G, (4.8)
where G ∼ N (0, 1) is a standard normal distribution.

The signal-to-noise ratio (SNR) can be dened by the ratio of mean to standard deviation of the measurement [Schroeder, 1999, Bushberg et al., 2003]

SN R i = λ i σ i = λ i √ λ i = √ λ i .
(4.9)

The above equation indicates that the SNR of the detected signal in each energy bin is proportional to the square root of the expected counting. Consequently, image noise in spectral CT or CR is inuenced by the number of photons falling inside each energy bin. Increasing the number of energy bins N could bring better resolution for material decomposition, but can also increase image noise. The trade-o between hiring more energy bins and augmenting SNR needs to be considered [START_REF] Leng | Noise reduction in spectral ct: Reducing dose and breaking the trade-o between image noise and energy bin selection[END_REF]. Another factor that directly inuences the number of photons arriving on the detector is the exposure factor (milliamperes×time). Therefore, in the present work, we investigated the inuence of these two factors on the performance of material determination methods, details will be presented in Chapter 5.

Decomposition quality assessment

The performance of material decomposition was evaluated by the average density, average error and the contrast-to-noise ratio (CNR) [START_REF] Leng | Maximizing iodine contrast-to-noise ratios in abdominal ct imaging through use of energy domain noise reduction and virtual monoenergetic dual-energy ct[END_REF].

Average density

For a region of interest (ROI) on the decomposed image, the average density ρROI can be dened by averaging all values computed after material decomposition in that region:

ρROI = ∑ J j=1 ρ j J , ( 4.10)
where ρ j denotes the calculated density of pixel j, and J the total number of pixels within that ROI.

Average error

Average error AE ROI (%) represents the mean value of absolute relative errors of the ROI, with respect to theoretical density t j :

AE ROI = 1 J J ∑ j=1 |ρ j -t j | t j .
(4.11)

Contrast-to-noise ratio

To calculate the CNR of the region of interest, we choose another region on the neighboring background area and denote it as BG. The mean density and the standard deviation (i.e.

image noise) of the BG region can be denoted respectively as ρ BG and SD BG . The CNR is then given by:

CNR = ρ ROI -ρ BG SD BG . ( 4.12)
As in our case, ρ ROI should be greater than ρ BG , we therefore dene the minimum and maximum CNR for the same ROI and BG as: Cardiovascular disease (CVD) produces immense health and economic burdens globally and atherosclerosis is a major contributor for CVD [START_REF] Benjamin | Heart disease and stroke statistics-2017 update: a report from the american heart association[END_REF]. Atherosclerosis is a disease in which the inside artery narrows due to the build up of plaque. Plaques inside the coronary arteries is a serious risk factor for adverse cardiovascular events. The ruptures of vulnerable plaques are the cause of about 70% of fatal acute myocardial infarctions and sudden coronary deaths [START_REF] Naghavi | From vulnerable plaque to vulnerable patient: a call for new denitions and risk assessment strategies: Part i[END_REF]. It is of great importance to detect the plaque and assess its vulnerability.

CNR min = ρ ROImin -ρ BGmax SD BG .
A high-risk plaque (HRP) usually demonstrates a large lipid-rich necrotic core that occupies a large proportion of the plaque, a thin brous cap with accumulated macrophages behind it, increased plaque inammation, neovascularization, luminal stenosis and so on [Moreno, 2010]. The commonly used invasive coronary imaging methods include intravascular ultrasound (IVUS), optical coherence tomography (OCT) and near infrared spectroscopy (NRIS) [START_REF] Sandfort | Noninvasive imaging of atherosclerotic plaque progression: status of coronary computed tomography angiography[END_REF]. Contrast-enhanced IVUS examination enables visualization of neovasculature of the atherosclerotic plaque and associated adventitial vasa vasorum [Feinstein, 2006]; the brous cap thickness can be assessed by OCT and the lipid blue, the lumen in green and the necrotic core in red. The larger the plaque area and the larger the necrotic core size, the higher is the likelihood of plaque vulnerability. Figure

retrieved from [Narula, 2009].

content of plaques can be assessed by NIRS [START_REF] Moreno | Detection of lipid pool, thin brous cap, and inammatory cells in human aortic atherosclerotic plaques by near-infrared spectroscopy[END_REF]. Noninvasive methods include magnetic resonance imaging (MRI) and radionuclide imaging, such as single photon emission computed tomography (SPECT) or positron emission tomography (PET).

MRI can demonstrate the coronary vessel lumen and characterizes the coronary vessel wall [START_REF] Macedo | Mri detects increased coronary wall thickness in asymptomatic individuals: The multi-ethnic study of atherosclerosis (mesa)[END_REF], Miao et al., 2009]; Technetium-99m-labelled ( 99m Tc) SPECT has the ability to detect the presence of inammation in carotid plaques [START_REF] Sun | Molecular imaging of plaques in coronary arteries with pet and spect[END_REF];

PET imaging can track plaque inammation using 18 F-uorodeoxy glucose (FDG) [START_REF] Rudd | Atherosclerosis inammation imaging with 18f-fdg pet: carotid, iliac, and femoral uptake reproducibility, quantication methods, and recommendations[END_REF]. All methods can detect the HRP to some extent with limited abilities and improved evaluation can be performed by a combination of multiple methods. However, it will be favorable to nd one imaging modality which can evaluate multiple characteristics of HRP at the same time and enable identication with better accuracy.

Computed tomography angiography (CTA) is an X ray computed tomography technique used to visualize arterial and venous vessels throughout the body. The patient is injected with contrast agent (usually iodinated contrast) and then the heart is scanned using a high speed CT scanner. There has been increasing interest in using CTA for atherosclerosis imaging. The multi-slice CT (MSCT) with conventional energy integrating detector has been widely implemented clinically and many research has been done to investigate its potential to detect coronary plaques. Preliminary studies show that MSCT can not only characterize lumen stenosis, but also quantify plaque volume and the extent of vessel wall re-modeling, and depict tissue density according to dierent X-ray attenuating properties of each structure [START_REF] Leber | Quantication of obstructive and nonobstructive coronary lesions by 64-slice computed tomography: a comparative study with quantitative coronary angiography and intravascular ultrasound[END_REF], Cademartiri et al., 2006, Voros et al., 2011]. Its performance for HRP detection has been validated by comparing with IVUS and OCT [START_REF] Voros | Prospective, head-to-head comparison of quantitative coronary angiography, quantitative computed tomography angiography, and intravascular ultrasound for the prediction of hemodynamic signicance in intermediate and severe lesions, using fractional ow reserve as reference standard (from the atlanta i and ii study)[END_REF], Nakazato et al., 2014]. However, MSCT characterizes composition of plaques depending on the absolute attenuation values, which are easily inuenced by multiple parameters: lumen attenuation, body mass index of the patient, and contrast to noise ratio of the images.

Dual-energy CT (DECT) can obtain additional spectral information compared with conventional energy integrating CT and has the potential to improve plaque imaging by enhancing identication of plaque compositions. Barreto et al. showed that DECT is able to characterize atherosclerotic plaques and distinguish calcied plaques from noncalcied plaques by comparison of attenuations at two tube voltages [START_REF] Barreto | Potential of dual-energy computed tomography to characterize atherosclerotic plaque: ex vivo assessment of human coronary arteries in comparison to histology[END_REF].

Tran et al. used DECT to identify the iodinated contrast and hydroxyapatite (HA, Ca 5 (PO 4 ) 3 (OH)), results show that DECT can detect and dierentiate between contrast medium and calcied tissues, but its accuracy is dependent on the CT density of tissues and limited when CT attenuation is low [START_REF] Tran | Dual-energy ct discrimination of iodine and calcium: experimental results and implications for lower extremity ct angiography[END_REF]. Bhavane et al. performed atherosclerosis imaging in mice with liposomal-iodine nanoparticle contrast agent to localize macrophage cells, and identied iodine and calcium concentration map based on DECT material decomposition [START_REF] Bhavane | Dual-energy computed tomography imaging of atherosclerotic plaques in a mouse model using a liposomal-iodine nanoparticle contrast agentclinical perspective[END_REF].

X-ray spectral CT obtains spectral information corresponding to dierent energy bins with single acquisition. It outperforms DECT with better material decomposition ability and shorter scan time, i.e. reduced radiation dose. Previous CTA studies have shown the limited ability to identify calcication of plaque and iodinated contrast agent which is used to assess stenotic narrowing. In the present study, we will investigate the decomposition ability of the methods (proj-LS, proj-LLS and proj-PR-LLS method) discussed in Chapter 3, and evaluate their performance in the identication of iodine, calcium and soft tissues through a computational human thorax phantom study.

Computational human thorax phantom

We use a computational human thorax phantom to simulate the in vivo situations. The geometry of the thorax phantom comes from VHP-Female v.2.2, which is a platformindependent full-body computational human model [START_REF] Noetscher | Vhp-female v3.0 fem/bem computational human phantom[END_REF]. This 3D CAD model contains up to 270 parts, 45 individual muscles, three body shells and so on, all extracted from open-source Visible Human Project -Female dataset of the National Library of Medicine. Figure 4.6 is an illustration of this model, where the geometries of dierent organs, bones, muscles, tissues, vessels can be clearly seen.

Our simulation needs only the thorax part of the full-body model and the blue line in , 1989] and earlier publications by White [START_REF] White | Average softtissue and bone models for use in radiation dosimetry[END_REF]. Based on these researches, we dened the components and densities of dierent tissues included in the thorax phantom shown in 4.7 (a). Table 4.1 is a list of this information. The simulated X-ray spectral CT system hired a 120kVp X-ray spectrum from tungsten target material without ltration, the target angle was 17 • and tube current was 100mA, as show in Figure 4.3. A linear cadmium telluride (CdTe) detector of 1000 pixels with pixel size of 1mm×1mm was used, and its detecting energy range was set from 30keV to 90keV with evenly distributed 6 energy bins. The absorption eciency of the detector was simulated for a 3 mm thick CdTe crystal. The thorax phantom was scanned with 1200 projections from 0 • to 360 • and the scan time was 0.01 s for each projection. The distances of source to the center of rotation and detector to the center of rotation were D SC = 700 mm and D DC = 300 mm respectively. Figure 4.8 illustrates the system geometry.

Results

In this application, we chose soft tissue, calcium and iodine as material basis during decomposition. Three projection domain material decomposition methods with dierent objective functions were applied:

1 ⃝ proj-LS method, corresponding to objective function in Equation (3.11).

2

⃝ proj-LLS method, corresponding to objective function in Equation (3.15).

3

⃝ proj-PR-LLS method, corresponding to objective function in Equation (3.16). Patch size was set to be 2 × 2 and relaxation parameter was chosen to be 1.

Figure 4.9 shows the decomposed images using the three methods. It is observed that all methods give the distribution of soft tissue, calcium and iodine:

• Skin, adipose, average soft tissue, cartilage, aorta, lung and heart tissues are visible only in the soft tissue image. Besides, part of the bones (ribs, vertebra and sternum) is also present in soft tissue image, because the components of bones contain not only calcium, but also other elements like hydrogen and oxygen.

• Calcium image well highlights bones and inserted calcium solutions

• Iodine is separated appropriately in the iodine-specic image.

We can observe obvious streak artifacts in the decomposed soft tissue image using proj-LS method, which is caused by individual acquisition pixels where few photons were detected and high noise resulted in visible measurement errors. This situation is much better when using the proj-LLS and proj-PR-LLS method. To make a quantitative analysis of the performance of the three methods, we calculated the average errors (AEs) and CNRs of the decomposed calcium and iodine solutions, as

shown in Figure 4.11 and Figure 4.12. ROIs were chosen as regions of calcium or iodine insert in the corresponding decomposed image (Figure 4.10 red circles). AE is calculated for each ROI by averaging the absolute relative errors using Equation (4.11). For the calculation of CNRs, we chose a region neighboring to the ROIs as the BG region (Figure 4.10 cyan circles), and determined CNRs according to Equation (4.14).

For the quantication of calcium, it is observed from Figure 4.11 (a) and Figure 4.12 (a)

that the proposed proj-PR-LLS method has the best quantication accuracy, i.e. lowest

AEs and best image quality, i.e. highest CNRs. For the determination of iodine (Figure 4.11 (b) and Figure 4.12 (b)), the proposed two methods have better performance than the commonly used proj-LS method. Comparing proj-LLS and proj-PR-LLS method, we can observe that proj-LLS method has lower quantication errors, while proj-PR-LLS method can generate iodine image with higher CNRs. By all the above analysis, we can conclude that our proposed proj-LLS and proj-PR-LLS method have better decomposition performance than the commonly used proj-LS method. Moreover, proj-PR-LLS method has the ability to enhance CNRs of decomposed images while obtaining better (calcium) or inferior but comparative (iodine) quantitative accuracy, when comparing with proj-LLS method. For visual benets, Figure 4.14 demonstrates the color overlay image using three decomposed material images shown in Figure 4.9 (c). The color overlay process was realized by setting the pixels with densities higher than 30 mg/cc in calcium image into red, the pixels with densities higher than 4 mg/cc in iodine image into green, and nally overlaying the two colors on the soft tissue image. Obvious streak artifacts can be observed in Figure 4.14, because there exists individual acquisition pixels where few photons were detected and high noise resulted in visible measurement errors after reconstruction, and the coloring process by setting a threshold made this artifact more obvious.

Conclusion and discussion

Motivated by the need to distinguish compositions of atherosclerotic plaque, we simulated spectral CT imaging process with a realistic computational human thorax phantom, and investigated the material decomposition performance of the three methods introduced in Chapter 3: proj-LS, proj-LLS and proj-PR-LLS method. Three images were obtained using either method: soft tissue, calcium and iodine image. Results showed that all methods can give quantitative distribution of three materials, among which soft tissue image contains skin, adipose, average soft tissue, cartilage, aorta, lung and heart tissues; calcium image well highlights bones and inserted calcium solutions; iodine image contains only the inserted iodine solutions. Quantitative analysis, including average errors and CNRs demonstrated that the proposed proj-LLS method and proj-PR-LLS method have better quantitative accuracy and better image quality than the commonly used proj-LS method.

Moreover, proj-PR-LLS method has the ability to enhance CNRs of decomposed images while obtaining better (calcium) or inferior but comparative (iodine) quantitative accuracy, when comparing with proj-LLS method.

There are potential improvements to our work. Radiation dose to patient is related to the number of projections. A lower number of projections can reduce patient dose but will also decrease image quality due to the sparser sampling of imaging data [Van Daatselaar et al., 2004]. The typical number of projections used by cone beam CT in applications of the breast, extremities, and head scans are 350-450 [START_REF] Zhao | Noise, sampling, and the number of projections in cone-beam ct with a at-panel detector[END_REF]. In the present study, 1200 projections were used in the simulated phantom scan for better image quality.

This number could lead to a prohibited patient dose and this should be considered in future studies.

Industrial context

Waste electrical and electronic equipment (WEEE) has been increasing rapidly due to the development of electronic industry. In European Union, the amount of WEEE generated in 2005 is 9 million tonnes and this number is supposed to grow to 12 million by 2020 [START_REF] Beigbeder | Study of the physico-chemical properties of recycled polymers from waste electrical and electronic equipment (weee) sorted by high resolution near infrared devices[END_REF], EuropeanCommission, 2018]. WEEE covers a wide variety of products such as lamps, hairdryers, computers, TV-sets, fridges and cell phones. This kind of waste contains various compositions among which polymers play an important part [START_REF] Cui | Mechanical recycling of waste electric and electronic equipment: a review[END_REF]. Acrylonitrile-butadiene-styrene (ABS) is a major component among all the polymers. The other polymers includes ABS/polycarbonate (PC) alloys, high-impact polystyrene (HIPS), polypropene (PP), poly(vinyl chloride) (PVC), and so on [START_REF] Barthes | Recycling of aged abs from real weee through abs/pc blends in the abs-rich compositions[END_REF]. Due to the presence of electronic power of EEE, the plastic materials should meet high re safety standards. However, it is not possible to realize with pure polymers, therefore ame retardant (FR) are added to change ammability of plastics and increase the re resistance [START_REF] Peeters | Closed loop recycling of plastics containing ame retardants[END_REF]. There are mainly four categories of FRs: (a) halogenated organic compounds (e.g., brominated aromatic compounds, chlorinated parans and alicyclic compounds), (b) phosphorus compounds (e.g., organophosphates, halophosphates, phosphine oxides and red phosphorus), (c) nitrogenbased compounds and (d) Inorganic salts [START_REF] Fromme | Brominated ame retardantsexposure and risk assessment for the general population[END_REF].

Recycling of plastics from WEEE is challenging because of the existence of FRs (especially the brominated and chlorinated FRs), which might result in serious environmental pollution. Therefore, the sorting technique is important for subsequent treatment and recycling of WEEE to avoid major environmental and health problems. Multiple techniques were applied for the plastics identication, such as: Near-infrared (NIR) analysis, Raman absorption spectrometry (RAS), Fourier-transform infrared (FTIR), sliding spark spectroscopy (SSSP), X-Ray uorescence (XRF), X-Ray computed radiography (CR), and dual Energy computed radiography (DECR) [START_REF] Peeters | Closed loop recycling of plastics containing ame retardants[END_REF].

CR technique measures the absorption of X-ray beam through objects. For materials with similar thickness, it is possible to separate them according to their attenuation nature to X-rays. Material with higher density and larger atomic number tends to absorbs more photons. This method can dierentiate metal and non-metal materials or lead and barium glasses. However, it is not possible to distinguish materials with similar attenuation properties. DECR provides extra spectral information for material decomposition. Research of using DECR for multi-material identication can be found in waste management [START_REF] Mesina | Automatic sorting of scrap metals with a combined electromagnetic and dual energy x-ray transmission sensor[END_REF], Montagner, 2012, Gundupalli et al., 2017]. For example, Duvillier et al. plotted the log-measurement using DECR, polyethylene and PVC can be separated from the curves [START_REF] Duvillier | Inline multi-material identication via dual energy radiographic measurements[END_REF]].

As we have introduced before, photon counting detector is able to obtain spectral information of several energy bins with single acquisition. It provided current CT technology with stronger material decomposition ability. However, the PCDs are not restricted to CT imaging, but can also be applied to CR. As we have mentioned in Chapter 2, the rst commercial photon-counting system was introduced for mammography [Åslund et al., 2007], which is an application of X-ray transmission imaging. Furthermore, the proposed patchwise regularized method in projection domain process directly the attenuated X-ray information, and can be applied for CR imaging. In the present work, we will investigate the material decomposition ability of proj-PR-LLS method for spectral CR imaging in application for ABS-FR plastic sorting. Furthermore, we also develop a classication based decomposition method for this specic application, named proj-PR-LLS-clas method. Its performance will be evaluated in comparison with proj-PR-LLS method.

ABS-FRs phantom

As we have said in the above section, ABS is a major component among all the polymers used by the electrical and electronic equipments, thus we selected it as the polymer material for investigation. We also chose three kinds of commonly used ame retardants, including the brominated (Br), chlorinated (Cl) and phosphorus (P) FRs. They are respectively tetrabromobisphenol A (TBBPA), dechlorane plus (DDC-CO) and resorcinol bis(diphenyl phosphate) (RDP). Three ABS-FR materials were obtained by mixing each FR with ABS at mass percentage of 15%. Information of these materials is summarized in Table 4.2. For example, ρ eff (Br) = ρ(mixture) × ω(Br) = 1060 × 8.82% = 93.5 mg/cm 3 ; ρ eff (ABS) = ρ(mixture) × ω(ABS) = 1060 × (1 -15%) = 901 mg/cm 3 .

We simulated a phantom composed of multiple cubes with height of 10 mm, width of 10 mm and dierent thicknesses. As shown in Figure 4.15, each column of cubes are of the same material (denoted on the top) and each row of cubes are of the same thickness (denoted on the left). Data taken from [START_REF] Hubbell | Tables of x-ray mass attenuation coecients and mass energy-absorption coecients (version 1.4)[END_REF]].

We will use projection domain decomposition method for this application, where the nal calculated values P α for material α are the line integral of material density ρ α (Equation 3.10). In this case, with known thickness d of each cube, theoretical values for P α at position (sx, sy) of the transmitted image can be given by P α (sx, sy) = d(sx, sy) × ρ α (sx, sy), (4.15) where ρ α equals the eective density calculated in Table 4.2. Table 4.3 lists the theoretical values of P α for material ABS, Br, Cl and P at various thicknesses.

System geometry

The simulated system uses the 100 kVp X-ray spectrum in Figure 4.3 with tube current of 15 mA. The distance from X-ray source to phantom center is 2000 mm, which is far enough to consider the X rays passing through the object as parallel. A 90×112 CdTe detector array with pixel size of 0.5 mm × 0.5 mm and thickness of 3 mm was simulated, and the distance of detector to phantom center was 3 mm. Six energy bins were set to be evenly distributed from 30 keV to 90 keV. Figure 4.17 is a schematic view of the simulated spectral CR system.

Planar array PCD

D DC =3 mm D SC =2000 mm X-ray source

Decomposition results of proj-PR-LLS method

For this application, we have simulated two sets of spectral CR acquisition data: with and without photon noise. Then we used proj-PR-LLS method to decompose the radiographic images into three basis images corresponding to ABS, Br and Cl. The reason for having not selected P as the basis is that P and Cl have too close atomic numbers and experiments showed that they can be barely separated. Thus we use Cl as basis material and expect that P will also be present in the Cl basis image. We set the patch size to be 2×2.

For simulated acquisition without noise, we set the relaxation parameter r to 0 since it is not necessary to enforce smoothness within patches when no noise exists, corresponding decomposition results are shown in When Poisson noise is added during the simulation, we need to reconsider the relaxation parameter r. r is determined as 10 4.4 in this situation according to the L-curve criterion. Some details of r selection using the L-curve method are introduced in Appendix section.

The decomposition results when r = 0 and r = 10 4.4 are shown in Figure 4.18 (b) and (c).

It can be seen that when r = 0, the ABS-TBBPA cubes and ABS-DDC-CO cubes appears in both Br basis image and Cl basis image, making them impossible to be distinguished.

Moreover, the decomposed images suers from heavy noise. When the relaxation parameter changes to r = 10 4.4 , the obtained images are less noisy than those in the former case, however, the separation of dierent materials is not enhanced.

Proj-PR-LLS-clas method

We have shown in the above section the decomposition results using proj-PR-LLS method.

Its performance is not satisfying, especially when Poisson noise is added during the simulation. To improve this situation, we develop another method based on proj-PR-LLS method.

As we have introduced, the objective of this application is to separate the ABS materials with and without ame retardant, and further identify the type of ame retardant.

In industrial applications, each plastic material contains only one type of FR and no overlapping exists between samples we perform CR inspection. Thus ideally, one pixel on the radiographic image of our phantom should either be identied as Br or Cl if it contains FR. Specically for this application, we develop a classication-based (proj-PR-LLS-clas) decomposition method. 4.3) of 3 basis materials and our calculations using two methods is illustrated. It is observed that the proj-PR-LLS-clas method has higher consistence with true values than the other one, especially for the Cl basis curve, where ABS-TBBPA cubes are mistaken as Cl contained FR by proj-PR-LLS method.

With simulated Poisson noise, the decomposed images of proj-PR-LLS-clas method also has higher material separation ability than those of proj-PR-LLS method. When r = 0, it can be observed from To conclude, the proj-PR-LLS-clas method has better material decomposition ability than proj-PR-LLS for this application, no matter when Poisson noise is or is not considered in the simulation. ABS and Br contained FR can be identied according to corresponding basis image. Furthermore, when the cube thickness increases to 2 mm and 4 mm, the decomposition results become much closer to those of noise-free condition, where ABS and the FRs containing Br, Cl and P can be identied at the same time. 4.3) if the decomposition is perfect, blue curves represent calculated values using method without classication and magenta curve using method with classication.

Conclusion and discussion

WEEE has been increasing rapidly due to the development of electronic industry. The sorting technique is important for subsequent treatment and recycling of WEEE to avoid major environmental and health problems. For this application, we have simulated the spectral CR imaging with a phantom containing multiple cubes of dierent ABS-FR materials with various thicknesses, and have investigates the decomposition feasibility of dierent ABS-4.6. SUMMARY FRs using the proj-PR-LLS method and proj-PR-LLS-clas method. Results show that the former method failed to identify the FR materials when Poisson noise was simulated for CR imaging. The latter method, in noise-free condition, could distinguish ABS and the FRs containing element Br, Cl and P at the same time, with quantication of each material agreeing well with the theoretical value. When noise was added in the simulation, the proj-PR-LLS-clas method can separate four materials at the same time only when the tube thickness is as large as 2 mm or 4 mm.

There are certain limitations in our work. Firstly, the two methods are not very robust to noise when the object thickness is too small, especially the proj-PR-LLS method.

Secondly, proj-PR-LLS-clas method has been proposed under the premise that each pixel on the radiographic image of object contains either one set of basis materials (ABS+Br) or the other (ABS+Cl), so that the material not belonging to the selected basis can be put 0. Proj-PR-LLS-clas method can not be applied for situation where this premise is not met. Thirdly, the decomposition of elements with similar atomic numbers (Cl and P), even in the noise-free condition of simulation, depends on their density dierence in the same image. Therefore, in case where the thickness information of objects is unknown, or the two materials overlap in the X-ray direction, the decomposed density will be largely inuenced, leading the materials to be indistinguishable.

Summary

In this chapter, we have presented the spectral CT/CR simulation method using the VXI software. The simulation is based on ray-tracing techniques together with the X-ray attenuation law. The users can dene parameters of the X-ray source, object and detector for their own applications. Furthermore, spectro mode of detector allows counting the arriving photons into corresponding pre-dened energy bins to simulate the spectral X-ray imaging process.

With VXI software, we have simulated spectral X-ray imaging with dierent phantoms, and investigated the material decomposition performance of the methods proposed in Chapter 3 for medical and industrial applications:

1. Medical application: iodine and calcium identication. To distinguish compositions of atherosclerotic plaque, we have simulated spectral CT imaging of a computational human thorax phantom with realistic noise and investigated the material decomposition performance of three methods: proj-LS, proj-LLS and proj-PR-LLS method.

Results show all methods can give quantitative distribution of iodine, calcium and soft tissues respectively, among which the proposed proj-LLS method and proj-PR-LLS method have better quantitative accuracy and better image quality than the commonly used proj-LS method.

3. Industrial application: ABS-FR plastic sorting. In addition to spectral X-ray CT, we have also investigated the material decomposition ability of the proposed proj-PR-LLS method for spectral CR imaging in the identication of dierent ABS-FR materials. Results show that the decomposition directly in the 3-matieral basis is not enough robust to noise. Thus we have introduced a classication-decomposition procedure (proj-PR-LLS-clas method) based on the application because we know a priori that one piece of plastic contain can not contain 2 FRs at the same time. So it is reasonable to rstly perform two decompositions in 2-materials basis and then choose the more likely one to be the nal result. This method has enhanced the separation between Br and Cl compared to proj-PR-LLS method. For simulation with realistic noise, when the thickness of phantom is favorable, ABS and the FRs that contain Br, Cl and P can be identied at the same time.

Page 81

Cette thèse est accessible à l'adresse : http://theses.insa-lyon.fr/publication/2018LYSEI056/these.pdf © [T. Su], [2018], INSA Lyon, tous droits réservés

Introduction

We have described in the previous chapter the material decomposition performance of X-ray spectral CT and spectral CR with photon counting detector. In the above applications, the spectral CT and CR imaging process were simulated using acquisition parameters commonly used in practical situations. For example, the number of energy bins was set to 6 considering that most existing PCDs (Table 2.1) have 2 to 8 bins. In medical imaging applications (atherosclerosis imaging), the X-ray tube current and exposure time were set to be 100 mA and 0.01 s/projection, resulting in an exposure factor (tube current and time product) of 1 mA•s per projection. What will be the performance of material decomposition if we vary these acquisition parameters during the simulation?

As indicated by Equation 4.9 and recalled here

SN R i = λ i σ i = λ i √ λ i = √ λ i , ( 5.1)
the signal-to-noise ratio of energy bin B i , dened by the ratio of mean (λ i ) to standard deviation (σ i ) of measurement, is proportional to the square root of the expected count- ing. The exposure factor has a direct inuence on the number of photons arriving at the detector, and consequently on the level of image noise. Therefore, we believe that the decomposition performance should be dierent if we change the exposure factor.

For a given X-ray spectrum, the width of each energy bin narrows when we increase the total number of energy bins, as illustrated in Figure 5.1. A narrow energy bin will limit the number of photons received by the detector in each bin and result in decreased SNR.

However, larger number of energy bins provides more energy dependent measurements and brings better resolution for material decomposition. Therefore, by increasing the number of energy bins, there exists a trade-o between decreasing SNR and increasing spectral information for material decomposition. In the present work, we will simulate X-ray spectral CT imaging of a poly(methyl methacrylate) (PMMA) phantom lled with iron (Fe), calcium (Ca) and potassium (K) [Mavrogeni, 2009,Ibrahim andBowman, 2014]. However, these existing techniques all have their own limitations. Biopsy is an invasive procedure that may cause bleeding from the incision site, pain and infections, also it is limited by sampling errors since iron is not uniformly distributed inside organs [START_REF] Yoon | Eectiveness and limitations of core needle biopsy in the diagnosis of thyroid nodules: review of current literature[END_REF].

Serum ferritin test is less invasive, but it can be easily inuenced by many factors such as inammation, infection, hepatic dysfunction, ascorbate deciency, hemolysis, and ineffective erythropoiesis [Crosby, 1976]. SQUID is particularly more accurate for liver iron content (LIC) determination but not for other organs [START_REF] Kolnagou | Uses and limitations of serum ferritin, magnetic resonance imaging t2 and t2* in the diagnosis of iron overload and in the ferrikinetics of normalization of the iron stores in thalassemia using the international committee on chelation deferiprone/deferoxamine combination protocol[END_REF]. MRI is a useful imaging tool for the diagnosis of iron deposition, but has very long scan time and inability to quantify iron concentration greater than 300 mol/g due to susceptibility and rapid signal decay artifacts [Ibrahim andBowman, 2014, Fischer et al., 2011]. Therefore, it is necessary to develop new techniques for tissue iron content determination.

In the past years, numerous technologies of DECT and spectral CT were developed, and their clinical applications were increasingly investigated [Fredenberg, 2017]. As we have mentioned, DECT has limited ability in material decomposition by producing two datasets corresponding to two dierent X-ray spectra of the same anatomic region [START_REF] Danad | New applications of cardiac computed tomography: dual-energy, spectral, and molecular ct imaging[END_REF]. The potential of using DECT to detect the distribution of iron in tissues was investigated [START_REF] Fischer | Quantication of liver iron content with ctadded value of dual-energy[END_REF], Hazirolan et al., 2008, Luo et al., 2015, Panta et al., 2015]. which depends on an empirical calibration with various concentrations of the chosen basis materials, and experimental validation was performed for the decomposition of calcium chloride (CaCl 2 ) and iron chloride (FeCl 3 ). The results showed that CaCl 2 and FeCl 3 can be separated while the performance is subject to energy threshold conguration of acquisition [START_REF] Li | Imagebased material decomposition with a general volume constraint for photon-counting ct[END_REF].

To summarize the literature survey, various methods exist for iron content evaluation, but with dierent limitations. At the same time, many material decomposition methods have been proposed and dierent experimental systems are under development for spectral CT. However, there is a lack of feasibility analysis and method comparison for the use of spectral CT for iron content determination. Therefore, the present study aims to investigate the ability of spectral CT to quantitatively separate iron from other elements like calcium and potassium, by using the proposed patchwise regularized decomposition methods in image domain (proj-PR-LLS method) and in projection domain (ima-PR-LS method). We simulated an 80mm-diameter PMMA phantom for scanning. Inside the phantom, 31 insert holes with diameter of 6mm were lled with dierent solutions, as shown in Figure 5.2. The rst two rows represent pure water rows that contain 7 inserts of pure water, and the third row correspond to the mixture row. The 3 bottom rows contain respectively iron, calcium and potassium soluted by water at dierent concentrations (15-200mg/cc). Specic density value for each insert is annotated inside each circle on the gure with unit of mg/cc. Symbol `#' marks the mixture insert of iron (Fe), calcium (Ca), potassium (K) and water. For example, the rst hole of the mixture row consists of 15mg/cc iron, 15mg/cc calcium and 15mg/cc potassium, the second one contains 30mg/cc iron, 30mg/cc calcium and 30mg/cc potassium, and so on. keV. Data taken from [START_REF] Hubbell | Tables of x-ray mass attenuation coecients and mass energy-absorption coecients (version 1.4)[END_REF].

Simulated phantom & system geometry

Figure 5.4 shows the system geometry of the simulated spectral CT system. The simulated system hired a linear CdTe detector of 700 pixels with pixel size of 0.4mm×0.4mm.

The absorption eciency of the detector was simulated for a 3 mm thick CdTe crystal.

The X-ray spectrum has tube voltage of 100kVp (spectrum shown in Figure 4.

Results

As mentioned in Chapter 3, the number of basis materials must be equal to or smaller than the number of energy bins. Therefore, for a complex object composed of multiple materials, all the components can not be put in the material basis. In this case, the attenuation caused by components having similar eective atom numbers Z eff are supposed to be represented by the same basis material. Based on this hypothesis, to comprehensively describe the whole object, we choose three materials as basis according to their atomic numbers: light material, median material and heavy material. In the case of our phantom, according to Figure 5.3, the light material corresponds to PMMA or water, median material Ca or K, heavy material Fe. To enhance the dierence between heavy (Fe) and median materials, we choose K as median material for both proj-PR-LLS and ima-PR-LS methods. However, dierent light materials are chosen for the two methods (PMMA for proj-PR-LLS and water for ima-PR-LS) because of their dierent decomposition domain. When we perform material decomposition in projection domain, the handled data are acquired sinograms which correspond to the transmitted number of photons of each ray. It should be noticed that when the X-ray passes through the phantom, PMMA, as a carrier of all inserts, contributes to most of the X-ray attenuation. Therefore, we select PMMA as the light material for proj-PR-LLS method. The situation is much dierent when we decompose in image domain, since we deal with the reconstructed images which corresponds to the attenuation coecients map. For each pixel of the reconstructed image, the attenuation coecient is caused by the components present in that pixel. Considering that water is a good solvent that is more likely to exist in the mixture than PMMA, we choose water as the light material for ima-PR-LS method.

Proj-PR-LLS method

As analyzed before, for proj-PR-LLS mehtod, the three materials chosen as basis material are PMMA, potassium and iron. Relaxation parameter r is set to be 1 and the size of each patch is 5×5. The results have shown that, iron can be identied in all conditions using proj-PR-LLS method. However, the performance changes a lot when the acquisition [START_REF]Problem statement and objective[END_REF][START_REF][END_REF], [START_REF]Problem statement and objective[END_REF]10), [START_REF]Problem statement and objective[END_REF]15), and so on. Y-axis represents dierent iron concentrations. Decomposition results with exposure factor of 0.025 mA•s are not presented in this plot for the sake of reading convenience. the CNR depends largely on the exposure factor and it benets from a higher exposure factor. For the same exposure factor, the CNR varies with the number of energy bins. We can observe that for lower exposure factor, decomposition performance tends to benet from a lower number of energy bins. For example, when the exposure factor is 0.025 mA•s, 6 energy bins gain the highest CNR. When exposure factor keeps increasing to 100 mA•s, the number of energy bins has little inuence on the CNR performance.

Ima-PR-LS method

For ima-PR-LS method, we choose water, potassium and iron as basis materials. Relaxation parameter r is set to be 1 and the size of each patch for regularization is 5×5. With various acquisition parameters, the simulation results show that with ima-PR-LS decomposition method, iron can be separated from calcium and potassium only when the number of energy bins is 60. We show in Figure 5. It is observed from Figure 5.10 (a) that the water basis image (left) that PMMA and water (rst 7 inserts) are present. Moreover, in regions of inserts containing dierent solutions where water is supposed to be separated, the density is higher than expected, i.e.

density of water 1000 mg/cc, especially for the row of mixture inserts. In iron basis image (middle), we can not only observe the two rows of inserts that really contain iron, but also the rows of inserts containing calcium and potassium which are supposed to be present in the potassium basis image, however, the potassium basis image (right) carries barely information. Figure 5.11 is the 1-D prole along the red line in Figure 5.10, theoretical values are given by the black curve. We can see extra peaks around the 500th and 600th pixel where calcium and potassium are contained. Therefore, iron is not separated from calcium and potassium in this case. Figure 5.12 is the 1-D prole along the blue and magenta lines in Figure 5.10. It can be observed that the measured iron density in the Fe/water inserts agrees well with the theoretical value while the number is much higher for iron in mixture inserts, caused by the fact that calcium and potassium are mistaken as iron. In conclusion, with the acquisition parameters of N =6 and F =1 mA•s, the proposed patchwise regularized method in image domain is able to separate iron from water, but fails to separate iron from calcium and potassium.

Conclusion and discussion

To investigate the inuence of exposure factor and number of energy bins on the performance of material decomposition, we have simulated spectral CT imaging of a PMMA phantom with Fe, Ca and K solutions at dierent concentrations. The simulations were done under various acquisition conditions by varing the number of energy bins (6, 10, 15, 20, 30 and 60) and exposure factors per projection (0.025, 0.1, 1, 10 and 100mA•s). Then we performed material decomposition using proj-PR-LLS method and ima-PR-LS method, which have been proposed in Chapter 3.

Page 98

Cette thèse est accessible à l'adresse : http://theses.insa-lyon.fr/publication/2018LYSEI056/these.pdf © [T. Su], [2018], INSA Lyon, tous droits réservés Proj-PR-LLS decomposition method enabled us to discriminate iron from calcium, potassium and water under all conditions. According to the CNR plot, for lower exposure factor, decomposition performance tends to benet from a lower number of energy bins.

This nding agrees with the theoretical analysis that low dose acquisition suers from more serious noise problem (Equation (5.1)), and that the more the energy bins, the more serious this problem is. We can employ the detector with an appropriate number of energy bins, depending on the acquisition parameter.

Ima-PR-LS decomposition method is much more sensitive to the number of energy bins than proj-PR-LLS method. Iron in the mixture was not separated from calcium and potassium when N = 6, 10, 15, 20 and 30. The only successful case appeared when N = 60.

The possible reason is that this method considers material decomposition as a linear inverse problem, where the eective mass attenuation coecients µ m (i) of each energy bin B i form its coecient matrix. However, µ m (i) was estimated by Equation (3.13) and this estimation has brought deviations to the true coecient matrix. The larger the number of energy bins, the more accurate the estimated coecient µ m (i) is and the less deviation it causes. With 60 energy bins, ima-PR-LS method outperforms the other approach, with higher CNRs.

In the case N = 60 and F = 100mA•s, the average errors of estimated iron content varied from 35.6% to 9.01% for the concentrations from 15mg/cc to 200 mg/cc.

We have mentioned in Chapter 3 that the proposed patchwise decomposition may cause pattern defect, which can be observed in the zoomed views in Figure 5.17. For image domain method, this defect exist in the form of patch eect, where gray dierence can be observed from one patch to another (Figure 5.17 (b)); for projection domain method, pattern defect exist in the form like ring artifact (Figure 5.17 Beam hardening is a well-known eect in X-ray CT scan, which expresses the phenomenon that the polychromatic X-ray beam becomes harder" when it passes through the object, due to the attenuation of lower energy photons and the presence of only higher energy photons. Some research showed that spectral CT imaging has the ability to correct beam hardening by projection domain material decomposition [START_REF] Roessl | K-edge imaging in x-ray computed tomography using multi-bin photon counting detectors[END_REF],

or by virtual monochromatic imaging technique [START_REF] Neuhaus | Improvement of image quality in unenhanced dual-layer ct of the head using virtual monoenergetic images compared with polyenergetic single-energy ct[END_REF], Yu et al., 2012]. Image domain decomposition works on reconstructed images of dierent energy bins, where beam hardening artifacts can be introduced by conventional reconstructions. Wider energy bin width (small number of energy bins) leads to heavier beam hardening. When the number of energy bins increases to 60 in our case, the beam hardening eect can be considered negligible (as the variation of the attenuation coecient in one keV is negligible).

Both proj-PR-LLS and ima-PR-LS methods for spectral CT are able to discriminate and quantify iron from calcium, potassium and water under certain conditions. Their performances vary with the acquisition parameters of spectral CT. One can choose to employ one of the two methods according to the acquisition constraints: if the number of energy bins is large (60 for example), we suggest ima-PR-LS method for iron determination in view of its higher CNR and shorter running time; otherwise, we suggest the proj-PR-LLS method.

Certain points of the present study could be improved. We have investigated the inuence of the number of energy bins on the decomposition performance, during which the bins were assumed to be evenly distributed within the detecting energy range. However, the width and placement of each energy bin can also inuence the performance of material Acquainted with the general principles of spectral CT, we set up the present work following the above-mentioned three aspects:

(1) Material decomposition methods (Chapter 3).

We have classied existing material decomposition methods in two categories: eectbased and material-based methods. Our work was around the latter one, where the linear attenuation coecient is described as the contributions of a basis of materials.

Meanwhile, according to the type of data on which the decomposition operates, we have investigated the projection domain method and the image domain method.

The commonly used decomposition is based on least square criterion, named proj-LS and ima-LS method. However, the inverse problem of material decomposition is usually ill-posed and the X-ray spectral CT measurements suer from Poisson photon counting noises. The standard LS criterion can lead to overtting to the noisy measurement data. We have proposed a least log-squares criterion for projection domain method to minimize the errors on linear attenuation coecient: proj-LLS method. Furthermore, to reduce the eect of noise and enforce smoothness, we have proposed to add a patchwise regularization term to penalize the sum of the square variations within each patch, giving rise to the proj-PR-LLS and ima-PR-LS methods.

(2) Validation and assessment via simulation (Chapter 4).

The performance of the proposed projection domain material decomposition methods has been validated through simulated spectral X-ray imaging. The simulation was realized with VXI software for medical and industrial applications:

-Medical application: iodine and calcium identication. To distinguish compositions of atherosclerotic plaque, we have simulated spectral CT imaging of a computational human thorax phantom with realistic noise. Extra calcium and iodine solutions with dierent concentrations were placed inside the phantom.

Three material decomposition methods have been investigated: proj-LS, proj-LLS and proj-PR-LLS methods. By decomposing the acquired data into three basis materials, iodine, calcium and soft tissues, all methods enable us to obtain the density distributions of the corresponding basis materials. However, the proposed proj-LLS method and proj-PR-LLS method exhibit better accuracy and better image quality than the commonly used proj-LS method in terms of their average errors and CNRs.

-Industrial application: ABS-FR plastic sorting. Although the material decomposition methods were demonstrated in case of spectral CT imaging, the projection domain decomposition methods can also be applied to PCD-based spectral CR imaging. For the latter, we have investigated the material decomposition ability of the proposed proj-PR-LLS method for spectral CR imaging in order to identify dierent ABS-FR materials (FRs: brominated FR, chlorinated FR and phosphorus FR). The results show that the decomposition directly in the 3-material basis (ABS+ Br+Cl) is not enough robust to noise. To improve the performance, we proposed a proj-PR-LLS-clas method which integrates classication into the decomposition procedure. The idea of classication is based on the fact that each plastic material contains only one type of FR. Therefore, it is reasonable to decompose rstly in 2-material basis (ABS+Br and ABS+Cl),

and then choose the one with smaller objective function values as the nal results. This method has enhanced the separation between Br and Cl compared to proj-PR-LLS method. For simulation with realistic noise, when the thickness of phantom is large enough, ABS and the FRs that contain Br, Cl and P can be identied at the same time.

(3) Inuence of dierent acquisition parameters (Chapter 5).

We have simulated spectral CT imaging of a PMMA phantom with Fe, Ca and K solutions at dierent concentrations to investigate the inuence of exposure factor and number of energy bins on the performance of material decomposition. The simulations were done under various acquisition conditions by varing the number of energy bins (N =6, 10, 15, 20, 30 and 60) and exposure factors per projection (F =0.025, 0.1, 1, 10 and 100mA•s). The performance of proj-PR-LLS method and ima-PR-LS has been studied. Proj-PR-LLS decomposition method enabled us to discriminate iron from calcium, potassium and water under all conditions. According to the CNR plot, for lower exposure factor, decomposition performance tends to benet from lower number of energy bins. Ima-PR-LS decomposition method is much more sensitive to the number of energy bins than proj-PR-LLS method. Iron in the mixture was not separated from calcium and potassium when N = 6, 10, 15, 20 and 30. The only successful case appeared when N = 60.

To conclude, both proj-PR-LLS and ima-PR-LS methods for spectral CT are able to discriminate and quantify iron from calcium, potassium and water under certain conditions. Their performances vary with the acquisition parameters of spectral CT.

PERSPECTIVES

One can choose to employ one or the other depending on acquisition constraints: if the number of energy bins is large (60 for example), we suggest ima-PR-LS method for iron determination in view of its higher CNR; otherwise, we suggest the proj-PR-LLS method.

Perspectives

We have presented in the previous parts the decomposition performance of the proposed methods, and the inuence of dierent acquisition parameters on the methods. Several potential improvements can be carried out in the future.

• On the simulation -The absorption eciency of the CdTe-based photon counting detector has been simulated in the present work. However, our simulation has assumed a perfect detector energy resolution so that all the photons can be correctly counted for the corresponding energy bin. In practice, the detector response function of monochromatic beam has a broad spectrum and the energy resolution is limited.

Therefore, the detector energy resolution should be considered in the future.

-Our simulation was based on the Beer Lambert law only, i.e. without scattered photons though VXI software allows simulation with rst order scattering at a cost of much longer simulation time. Therefore, taking scattering into consideration for the simulation will be a step forward to see if it has an inuence on the decomposition performance.

-State of the art PCDs have defects including pulse pileup, charge sharing, K-escapes, Compton scattering and so on, which were not simulated in our study and can be improved in the future as published by others. In the work of [START_REF] Taguchi | An analytical model of the eects of pulse pileup on the energy spectrum recorded by energy resolved photon counting x-ray detectors[END_REF], Alvarez, 2014], pulse pileup eect has been modeled by analytical formalism. Furthermore, in the work published by CEA-LETI [START_REF] Montémont | Simulation and design of orthogonal capacitive strip cdznte detectors[END_REF], a model of the semiconductor detector has been developed which included other defects such as charge sharing, K-escapes and Compton scattering in the detector response matrix.

• On the decomposition methods -We have used FBP reconstruction throughout the present work for sake of computation cost. However, iterative reconstruction has potential to reduce noise and improve image quality. It would be interesting to replace FBP reconstruction with iterative methods if computational condition is allowed.

-For the proposed patchwise regularized methods (proj-PR-LLS and ima-PR-LS), several parameters need to be chosen: patch size (np r × np c ) and relaxation parameter (r). We have used single patch size for the whole image of each application, further study can be revolved around a zone-wise patch size.

For example, we take the low energy bin image and apply a segmentation process to detect the homogeneous zones and detail zones, then we dene large patch sizes for the former and small patch size for the latter. As for r, it has been chosen heuristically through visual inspection, except for the application of plastic sorting, where we have applied L-curve method to a certain region on the radiographic images, and used the approximated r for all patches. However, r should change for dierent patches since the noise level varies. Therefore, further improvement could be done concerning the selection of patchwise relaxation parameter.

• Experimental validation

The evaluation of the proposed methods in this thesis was based on simulated spectral X-ray imaging, which is dierent from real experimental systems. A Philips spectral CT prototype has been installed at imaging center CERMEP (Centre d'Etude et de Recherche Multimodal Et Pluridisciplinaire) in Lyon. It can detect and quantify a variety of atoms including iodine, gadolinium, gold and bismuth using K-edge technique, and also enables to enhance the contrast-to-noise ratio of images by virtual monochromatic imaging technique. However, the raw data are not available now for processing. Further experimental validation on real spectral CT data will be considered in the future.

Author's publications

Journal papers:

Appendix

Determination of r with L-curve method (section 4.5)

The L-curve is a log-log plot of the norm of a regularized solution versus the norm of the corresponding residual, where the corner of the curve corresponds to the selected relaxation parameter [Hansen, 1999]. For proj-PR-LLS decomposition method, the solution norm and residual norm correspond to the objective function value and R(P C α ) value in Equation 3.16. We choose a region of interest on the radiographic image for the determination of r, which is shown in The ROI contains 20×20 pixels. Known that the patch size has been set to 2×2, there are 10 ×10=100 patches within the ROI. For a given r, we dene the solution norm and residual norm as the average values calculated for these 100 patches. We plotted the Lcurve with r varying from 10 -10 to 10 10 , as shown in La tomographie (CT) aux rayons X joue un rôle important dans l'imagerie non invasive depuis son introduction. Au cours des dernières années, de nombreuses avancées technologiques en tomographie par rayons X ont été observées, notamment la CT spectrale, qui utilise un détecteur à comptage de photons (PCD) pour discriminer les photons transmis correspondant à des bandes d'énergie sélectionnées afin d'obtenir une information spectrale.

La CT spectrale permet de surmonter de nombreuses limitations des techniques précédentes et ouvre de nombreuses applications nouvelles, parmi lesquelles la décomposition quantitative des matériaux est le sujet le plus étudié. Un certain nombre de méthodes de décomposition des matériaux ont été rapportées et différents systèmes expérimentaux sont en cours de développement pour la CT spectrale. Selon le type de données sur lequel l'étape de décomposition fonctionne, nous avons les méthodes du domaine des projections (décomposition avant reconstruction) et les méthodes du domaine de l'image reconstruite (décomposition après reconstruction).

La décomposition couramment utilisée est basée sur le critère des moindres carrés, nommée proj-LS et méthode ima-LS. Cependant, le problème inverse de la décomposition du matériau est généralement mal posé et les mesures du CT spectral aux rayons X souffrent de bruits de comptage de photons de Poisson. Le critère des moindres carrés peut conduire à un surajustement des données de mesure bruitées. Dans le présent travail, nous avons proposé un critère de moindre log-carré pour la méthode du domaine de projection afin de minimiser les erreurs sur le coefficient d'atténuation linéaire: méthode proj-LLS. De plus, pour réduire l'effet du bruit et lisser les images, nous avons proposé d'ajouter un terme de régularisation par patch pour pénaliser la somme des variations au carré dans chaque zone pour les décompositions des deux domaines, nommées proj-PR-LLS et ima -PR-LS méthode.

Les performances des différentes méthodes ont été évaluées par des études de simulation avec des fantômes spécifiques pour différentes applications: (1) Application médicale: identification de l'iode et du calcium. Les résultats de la décomposition des méthodes proposées montrent que le calcium et l'iode peuvent être bien séparés et quantifiés par rapport aux tissus mous.

(2) Application industrielle: tri des plastiques avec ou sans retardateur de flamme. Les résultats montrent que 3 types de matériaux ABS avec différents retardateurs de flamme peuvent être séparés lorsque l'épaisseur de l'échantillon est favorable.

Enfin, nous avons simulé l'imagerie par CT spectrale avec un fantôme de PMMA rempli de solutions de Fe, Ca et K.

Figure 1 :

 1 Figure 1: Schéma de l'imagerie CT: l'image en coupe transversale est reconstruite à partir de données de projection mesurées à partir d'angles multiples..

Figure 2 :

 2 Figure 2: Illustration de l'acquisition de données à partir de diérents systèmes CT: système CT classique où les détecteurs intégrateurs acquièrent un seul sinogramme (a); DECT utilisant une technique de commutation rapide des kilovolts ou une technique à double source (b), ou une technique de détection à double couche (c) pour obtenir deux ensembles de données de projection; Le système de CT spectrale (d) avec un détecteur de comptage de photons qui acquiert plusieurs ensembles de données de bandes d'énergie sélectionnées.

  la validation et l'évaluation des méthodes proposées par simulation d'applications médicale et industrielle; (3) l'étude de l'inuence de diérents paramètres d'acquisition, c'est-à-dire le facteur d'exposition et le nombre de bandes d'énergie, sur les performances de décomposition. La thèse s'organise comme suit: Dans le Chapitre 2, intitulé General principles of spectral CT, les principes des interactions des rayons X avec la matière sont introduits, y compris l'absorption photoélectrique, la diusion Compton et la production de paires; Les principes de la CT par rayons X sont présentés, en ce qui concerne la CT conventionnelle, la CT à double énergie et la CT spectrale. Dans le Chapitre 3, intitulé Material decomposition methods, une revue des méthodes de décomposition en matériaux existantes pour la CT spectrale est présentée. Sur la base des méthodes existantes, nous avons proposé des méthodes de décomposition régularisées en zones (patch) dans le domaine des projections et dans le domaine de l'image reconstruite. Dans le Chapitre 4, intitulé Material decomposition results for dierent applications: phantoms study, la procédure de simulation d'imagerie spectrale aux rayons X est détaillée et les résultats de la décomposition par les méthodes proposées est démontrée par une étude Dans le Chapitre 5, intitulé Exploring the inuence of dierent acquisition parameters, plusieurs séries d'acquisitions de CT spectrale ont été simulées avec diérents paramètres et l'inuence de ces paramètres sur la performance de décomposition des méthodes proposées est étudiée. Dans le Chapitre 6, intitulé Conclusion and perspectives, un bref résumé des principaux résultats et perspectives d'avenir est présenté.

Figure 3 :

 3 Figure 3: Illustration de la décomposition du domaine de projection et d'image. L'exemple ici suppose que le détecteur CT spectral a 6 bandes d'énergie et donc 6 sinogrammes sont obtenus pour une acquisition. Le premier procédé décompose d'abord les données d'acquisition en M sinogrammes cibles (ici M = 4 matériaux) en fonction de la localisation, puis eectue une reconstruction sur chaque sinogramme obtenu pour obtenir l'image spatiale nale. La deuxième méthode reconstruit d'abord les images spatiales à partir des données de projection qui correspondent à certaines bandes d'énergie, puis eectue une décomposition en matériaux sur ces données dans le domaine image (i.e. les images spatiales reconstruites).

Figure 4 :

 4 Figure 4: Les méthodes de décomposition en matériaux discutées et proposées (en bleu) dans ce travail de thèse. Le terme CR fait référence à la radiographie spectrale.

Figure 5 :

 5 Figure 5: Coupe transversale du fantôme du thorax avec illustration de diérents matériaux inclus (a) et de la région du c÷ur grossie qui contient des inserts de calcium et d'iode avec des concentrations diérentes (b).

Figure 6 :

 6 Figure 6: Schéma d'acquisition du CT spectral avec fantôme humain.

Figure 7 :Figure 8 :Figure 9 :

 789 Figure 7: Résultats de la décomposition en matériaux avec la méthode proj-LS (a), la méthode proj-LLS (b) et la méthode proj-PR-LLS (c). Chaque colonne représente la même base de matériaux, i.e. le tissu mou pour la première colonne, le calcium pour la seconde et l'iode pour la troisième.

Figure 10 :Figure 11 : 2 × 2 ,

 101122 Figure 10: Illustration du fantôme ABS-RF utilisé pour l'imagerie CR spectrale.

Figure 12 :

 12 Figure 12: Les résultats de décomposition dans diérentes conditions en utilisant la méthode proj-PR-LLS. (a): simulation de CR sans bruit, r=0; (b): simulation de CR avec bruit de Poisson, r=0; (c): simulation de CR avec bruit de Poisson, r = 10 4.4 . La première colonne représente la base de l'ABS, la deuxième la base de Br et la troisième colonne la base de Cl. L'échelle de couleur représente le facteur d'atténuation (c'est-à-dire µ×epaisseur).

Figure 13 :

 13 Figure 13: Organigramme de la méthode de décomposition avec classication pour un patch. fval1 et fval2 représentent les valeurs des fonction de coût obtenues à la n de la décomposition 1 et de la décomposition 2.

Figure 14 :Figure 15 :Figure 16 :

 141516 Figure 14: Les résultats de décomposition dans diérentes conditions en utilisant la méthode proj-PR-LLS-clas. (a): simulation de CR sans bruit, r = 0; (b): simulation de CR avec bruit de Poisson, r = 0; (c): simulation de CR avec bruit de Poisson, r = 10 4.4 . La première colonne représente la base de l'ABS, la deuxième la base de Br et la troisième colonne la base de Cl.

Figure 17 :Figure 18 :

 1718 Figure 17: Illustration des composants du fantôme (a) pour mieux comprendre les résultats de la décomposition. Densité décomposée en unité de mg/cc en utilisant la méthode proj-PR-LLS avec paramètres d'acquisition: 6 bandes d'énergie, 1 mA•s par projection pour (b), et 60 bandes d'énergie, 10 mA•s par projection pour (c), où la colonne de gauche représente la base de PMMA, la colonne du milieu le fer , et la colonne de droite le potassium.

Figure 20

 20 Figure 20 une décomposition qui a échouée (a) N = 6, et F = 1 mA•s, et une décomposition réussie (b) lorsque N = 60, et F = 10 mA•s. On observe à partir de la Figure 20 (a) que dans l'image à base de fer (milieu), on peut non seulement observer les deux rangées d'inserts qui contiennent réellement du fer, mais aussi les rangées d'inserts contenant du calcium et du potassium qui ne devraient être présents que dans l'image de base de potassium, alors que l'image de base de potassium (à droite) porte à peine l'information. La Figure 20 (b) montre les résultats de la décomposition lorsque N = 60, and F = 10 mA•s. Une diérence évidente peut être observée par rapport à (a). Seulement les deux rangées d'inserts existent dans l'image de base de fer (au milieu), représentant le fer dans le mélange et dans l'eau; tandis que dans l'image du potassium (à droite), les inserts contenant du calcium et du potassium sont bien mis en évidence.

Figure 20 :Figure 21 :Figure 22 :

 202122 Figure 20: Densité décomposée en unité de mg / cc de base d'eau (colonne de gauche), base de fer (colonne du milieu) et base de potassium (colonne de droite) avec la méthode ima-PR-LS. Paramètres d'acquisition des données: 6 bandes d'énergie, 1 mA•s par projection (a), et 60 bandes d'énergie, 10 mA•s par projection (b).

2 . 7

 27 Dominant regions of three kinds of X-ray interactions with matter. Figure retrieved from website [Nuclear, 2018]. . . . . . . . . . . . . . . . . . . . . .

2. 10

 10 Spectrum of a tungsten X-ray source with tube voltage of 150 kV. The Bremsstrahlung emission, characteristic X-rays and the ltering of the lowerenergy part are clearly shown in the spectrum. Figure retrieved from website [Cattin, 2016] . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.11 Schematic description of a scintillator-based detector. . . . . . . . . . . . . . 2.12 Mechanics of the room-temperature semiconductor detector. . . . . . . . . . 2.13 Parallel projections P θ (t) of object f (x, y) for angle θ. . . . . . . . . . . . . 2.14 The Fourier slice theory: the Fourier transform of a set of projections along angle θ equals the Fourier transform of the object along the line AB. . . . . 2.15 The bandlimited lter response (a), and the corresponding impulse response (b). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.16 Fan beam projections where the detectors are positioned with equal space along a straight line. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.17 Parameters used in the derivation of the reconstruction method for eauispaced fan beam CT imaging. . . . . . . . . . . . . . . . . . . . . . . . . . . 2.18 Data acquisition from dierent CT systems: conventional CT system where integrating detectors acquire one single sinogram (a); DECT using rapid kilovoltage switching technique or dual-source technique (b) or dual-layer detector technique (c) acquiring two set of projection data; spectral CT system (d) with photon counting detectors acquire several datasets of selected energy bins. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.19 Principle of the photon counting detector. . . . . . . . . . . . . . . . . . . . 3.1 Energy dependence of X-ray attenuation coecient for iron (a) and iodine (b). K-edge of iron appears at 7.11 keV; iodine K-edge at 33.17 keV. Ledges of iodine : L-I 5.19 keV, L-II 4.85 keV and L-III 4.56 keV. Data taken from [Hubbell and Seltzer, 2004]. . . . . . . . . . . . . . . . . . . . . . . . . 3.2 Illustration of projection domain decomposition and image domain decomposition. The example here assumes that the spectral CT has 6 energy bins and 6 sinograms are obtained for one acquisition. The former method rstly decomposes acquisition data into M location-dependent target sinograms (M = 4) in this case, and then perform a reconstruction on each obtained sinogram to get the nal spatial image. Image domain method rstly reconstructs the spatial images from their respective projection data corresponding to a given energy bin, and then performs material decomposition on these image domain data (i.e. the reconstructed spatial images). . 3.3 Illustration of the proj-PR-LLS method for patch C. Patch C is noted by rectangular regions in the acquisition sinograms and the decomposed images. The objective function is calculated by summing up the squared residuals from acquisition data and the regularization term from the target decomposed data. Decomposed images in patch C is obtained by minimizing the objective function. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.4 Illustration of the image domain decomposition process based on a patchwise regularized least squares criterion for patch C. Patch C is denoted by the rectangular regions in the reconstructed or decomposed images. The objective function is calculated by summing up the squared residuals and the regularization term from the target decomposed data. Decomposed images in patch C is obtained by minimizing the objective function. . . . . . . . . .

4 . 2

 42 Illustration of the simulated imaging chain. Every box in this gure corresponds to a section of the simulation software, in which several parameters are adjustable. Figure is retrieved from[START_REF] Duvauchelle | A computer code to simulate x-ray imaging techniques[END_REF]. . . . . . . 54 4.3 The X-ray spectra are computed based on Birch & Marshall model for tungsten target material without ltration, the target angles were 17 • and the tube voltage were set to be 100kVp and 120kVp. . . . . . . . . . . . . . . . . 55 4.4 Scheme of sample denition process: geometry of each part of the sample is designed separately, then assembled together to obtain the sample geometry, nally we can dene the material for each part of sample to get the completed sample. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55 4.5 Histopathological characteristics of a ruptured plaque (A, B) and the plaque which is prone to disruption (C-D). In images B and D, the plaque area has been colored blue, the lumen in green and the necrotic core in red. The larger the plaque area and the larger the necrotic core size, the higher is the likelihood of plaque vulnerability. Figure retrieved from [Narula, 2009]. . . 59 4.6 Illustration of VHP-Female v2.2 model. Geometries of dierent body parts have been described in this model, including the organs, bones, muscles, tissues, vessels and so on. The blue line in this gure denotes the slice that is scanned during our later simulation. . . . . . . . . . . . . . . . . . . . . . 61 4.7 Cross section of the thorax phantom with illustration of dierent material types included (a) and the magnied heart region that contains extra calcium and iodine inserts with dierent concentrations (b). . . . . . . . . . . . 62 4.8 Scheme of spectral CT acquisition with computational human phantom. . . 62 4.9 Material decomposition results of proj-LS method (a), proj-LLS method (b) and proj-PR-LLS method (c). Each column represents the same material basis, i.e. soft tissue for the rst column, calcium for the second and iodine for the third. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64 4.10 Illustration of the ROI and BG regions chosen for quantitative analysis of calcium (a) and iodine (b). The red circles represent ROI and the cyan circles BG regions. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65 4.11 Average errors of decomposed densities of calcium (a) and iodine (b) using proj-LS, proj-LLS and proj-PR-LLS methods. . . . . . . . . . . . . . . . . 66 4.12 CNRs of decomposed calcium and iodine image at dierent densities using proj-LS, proj-LLS and proj-PR-LLS method. The error bars represent the variation of CNR from CNR min to CNR max . . . . . . . . . . . . . . . . . . . 67 4.13 Comparison of the decomposed average density of calcium (a) and iodine (b) with true values using proj-PR-LLS method. The error bars indicate the standard deviation of decomposed density. . . . . . . . . . . . . . . . . 68 4.14 Color overlay image of soft tissues(gray), calcium (red) and iodine (green). . 69 4.15 Illustration of the ABS-FR phantom used for spectral CR imaging. Materials of these cubes are introduced in Table 4.2. . . . . . . . . . . . . . . . . 72 4.16 Mass attenuation coecients (µ m ) of components (ABS, Br, Cl and P) con- tained in this phantom within the detecting energy range from 30 keV to 90 keV. µ m of ABS is calculated according to formula µ m = ∑ α ω α µ mα , with ω α and µ mα representing the weight fraction and mass attenuation coecient of each element α of the compound. Data taken from [Hubbell and Seltzer, 2004]. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.17 Scheme of spectral CR system geometry. This gure is only a simplied schematic view of the simulation system where the dimensional information (D SC , D DC , size of phantom and detector) is not in realistic scale. . . . . . 4.18 Decomposition results under dierent conditions using proj-PR-LLS method. (a): simulation of CR without noise, r = 0; (b): simulation of CR with Poisson noise, r = 0; (c): simulation of CR with Poisson noise, r = 10 4.4 . The rst column represents ABS basis, second column Br basis, and third column Cl basis image. The dash lines on images of (a) mark the position where 1-D proles are plotted in later analysis. . . . . . . . . . . . . . . . . 4.19 Flowchart of the decomposition method with classication for one patch. fval1 and fval2 represent the cost function values obtained at the end of decomposition 1 and decomposition 2. . . . . . . . . . . . . . . . . . . . . . 4.20 Decomposition results under dierent conditions using proj-PR-LLS-clas method. (a): simulation of CR without noise, r = 0; (b): simulation of CR with Poisson noise, r = 0; (c): simulation of CR with Poisson noise, r = 10 4.4 . The rst column represents ABS basis, second column Br basis, and third column Cl basis image. The dash lines on images of (a) mark the position where 1-D proles are plotted in later analysis. . . . . . . . . . . . 4.21 Performance comparison of proj-PR-LLS method and proj-PR-LLS-clas method in noise-free condition: 1-D proles along the dash lines in Figure 4.18 and 4.20. ABS basis (top), Br basis (middle) and Cl basis (bottom). Black curves represent the theoretical density integrals of basis materials (P α for thickness of 1 mm in Table 4.3) if the decomposition is perfect, blue curves represent calculated values using method without classication and magenta curve using method with classication. . . . . . . . . . . . . . . . . . . . . 5.1 Illustration of spectral CT acquisition with 4 (a) and 6 (b) energy bins for the same X-ray spectrum. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5.2 Phantom illustration. Each row contains the same material or materials soluted by water at dierent concentrations, types of material are indicated on the left hand and concentrations are annotated inside each insert on the gure with unit of mg/cc. Symbol `#' marks the mixture inserts which combine iron, calcium, potassium and water together. . . . . . . . . . . . .

5 . 4

 54 Scheme of spectral CT acquisition. . . . . . . . . . . . . . . . . . . . . . . . 5.5 Ground truth image of iron density distribution (a) and zoomed view of the red rectangle region (b). The circular regions of interest are placed over the 100mg/cc mixture solution and its neighboring background. . . . . . . . . . 5.6 Components illustration of phantom (a) to better understand the decomposition results. Decomposed density in unit of mg/cc using the proj-PR-LLS method with acquisition parameters: 6 energy bins, 1 mA•s per projection for (b), and 60 energy bins, 10 mA•s per projection for (c), where the left column represents PMMA basis, middle column iron basis and right column potassium basis. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 90 5.7 Comparison of the decomposed average density with true values by proj-PR-LLS method: iron (a) and potassium (b) when acquisition parameters are 6 energy bins, 1 mA•s per projection; iron (c) and potassium (d) when acquisition parameters are 60 energy bins, 10 mA•s per projection. The error bars indicate the standard deviation of decomposed density. . . . . . . 91 5.8 Average errors (AE) of decomposed iron (in the mixture) concentrations using proj-PR-LLS method. X-axis represents 4 groups of acquisitions corresponding to 4 dierent exposure factors F , in each of which are indicated 6 dierent numbers of energy bins N : for example, for the top right group,(F, N ) = (0.1, 6), (0.1, 10), (0.1, 15), (0.1, 20), (0.1, 30), (0.1, 60), (1, 6), (1, 10), (1, 15), and so on. Y-axis represents dierent iron concentrations. Decomposition results with exposure factor of 0.025 mA•s are not presented in this plot for the sake of reading convenience. . . . . . . . . . . 92 5.9 CNRs of decomposed iron images for 100 mg/cc of iron with dierent acquisitions using proj-PR-LLS method. Each group of acquisitions has the same number of energy bins but dierent exposure factors. . . . . . . . . . 92 5.10 Decomposed density in unit of mg/cc of water basis (left column), iron basis (middle column) and potassium basis (right column) by ima-PR-LS method. Data acquisition parameters: 6 energy bins, 1 mA•s per projection (a), and 60 energy bins, 10 mA•s per projection (b). The red, blue and magenta lines on the iron basis image mark the position where 1-D proles are plotted in later analysis. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 93 5.11 1-D proles along the red line in Figure 5.10 (a): black curve represents theoretical values, red curve represents measured values. . . . . . . . . . . 94 5.12 1-D proles along the blue and magenta lines in Figure 5.10 (a): black curve represents theoretical values, magenta curve iron/water inserts and blue curve mixture inserts. . . . . . . . . . . . . . . . . . . . . . . . . . . . 94 5.13 Relationship between calculated and true densities in the case of iron (a) and potassium (b) with true values by ima-PR-LS method. Data acquisition parameters: 60 energy bins, 10 mA•s per projection. . . . . . . . . . . . . . 96 5.14 Average errors comparison of decomposed iron images obtained using proj-PR-LLS method (black circles) and ima-PR-LS method (red diamonds). All acquisitions are with 60 energy bins but dierent exposure factors. . . . . . 97 5.15 CNRs comparison of decomposed iron images for 100 mg/cc of iron obtained using proj-PR-LLS method and ima-PR-LS method. All acquisitions are with 60 energy bins but dierent exposure factors. . . . . . . . . . . . . . . 97 5.16 Computation time of proj-PR-LLS and ima-PR-LS methods for dierent acquisitions: exposure factor F from 0.025 to 100 mA•s, energy bins N from 6 to 60 for proj-PR-LLS method and N = 60 for ima-PR-LS method. . . . 98 tained by proj-PR-LLS (a) and ima-PR-LS (b) methods to show the pattern defects caused by patchwise regularization. . . . . . . . . . . . . . . . A.1 Illustration of ROI for the determination of r. . . . . . . . . . . . . . . . . . A.2 The L-curve with r varying from 10 -10 to 10 10 (a), and from 10 4 to 10 5 (b).

4 . 1

 41 Elemental compositions and mass densities of body tissues included in the thorax phantom. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .4.2Components of three ABS-FR materials used for the phantom. . . . . . . .4.3Theoretical values of P α (unit: mg/cm 2 ) for each material. Material ABS has two columns of values, the rst is ABS pure, which corresponds to the material of the 4th column of cubes in the phantom shown in Figure 4.15, and ABS in mixture corresponds to the ABS components in the other cubes of Figure 4.15. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . X-ray computed tomography (CT) is rst introduced in the 1970s by Cormack and Hounseld, for which they were awarded the Nobel Prize for Physiology or Medicine in 1979. This invention showed the possibility to invasively give high quality cross-sectional images of object with multiple-angle measurements and more or less complex computations. Afterwards, the advancement of X-ray CT has greatly beneted from developed reconstruction algorithms, advanced computer technology and technological improvements of detectors and sources. It now has wide applications in medical, industrial, material, aerospace and geology elds.

Figure 1 . 1 :

 11 Figure 1.1: Contributions of the present work.

Figure 2 . 1 :

 21 Figure 2.1: The X-ray picture of the hand of Röntgen's wife, retrieved from his original publication [Röntgen, 1896].

  Figure 2.2 shows the wavelength, frequency and energy of electromagnetic spectrum. It can be observed that X-rays have high energy among the spectrum.

Figure 2

 2 Figure 2.2: The electromagnetic spectrum, presented as a function of wavelength, fre-

Figure 2 . 3 :

 23 Figure 2.3: During the photoelectric event, the incident photon is annihilated. Part of its energy is consumed by ionization of atom and the rest is given to ejected photoelectron as its kinetic energy. The transfer of electron from outer shell to inner shell to ll in the vacancy created by photoelectron will generate either a. X-ray uorescence, or b. Auger electron. Figure is inspired by[START_REF] Barrett | Radiological imaging: the theory of image formation, detection, and processing[END_REF].

Figure 2 . 4 :

 24 Figure 2.4: The mechanism of Compton scattering: incident photon is scattered, its energy and direction change and a recoiling electron is ejected. Figure is inspired by[START_REF] Barrett | Radiological imaging: the theory of image formation, detection, and processing[END_REF].

Figure 2 . 5 :

 25 Figure 2.5: The mechanism of pair production: incident photon is annihilated, an electronpositron pair is generated. Figure is inspired by[START_REF] Barrett | Radiological imaging: the theory of image formation, detection, and processing[END_REF].

Figure 2 . 6 :

 26 Photon energy (MeV)

  2.7. It is observed that pair production dominates only when photon energy is extremely high, above approximately 10 Mev. For lower energy which are commonly used in diagnostic CT, photoelectric eect and Compton scattering and Compton scattering has major inuence on low-Z materials.

Figure 2 .

 2 Figure 2.7: Dominant regions of three kinds of X-ray interactions with matter. Figure retrieved from website [Nuclear, 2018].

Figure 2 .Figure 2

 22 Figure 2.8 shows the mechanics of CT imaging, the basic components are X-ray source, object and detector.

Figure 2

 2 Figure 2.9: The EMI Mark I scanner (a), and the cross sectional image of brain (b).

Figure 2 .

 2 Figure 2.10: Spectrum of a tungsten X-ray source with tube voltage of 150 kV. The Bremsstrahlung emission, characteristic X-rays and the ltering of the lower-energy part are clearly shown in the spectrum. Figure retrieved from website[Cattin, 2016]

Figure 2 .Figure 2 .

 22 photodiode layer.

Figure 2 .

 2 Figure 2.13: Parallel projections P θ (t) of object f (x, y) for angle θ.

Figure 2 Figure 2 .

 22 Figure 2.14: The Fourier slice theory: the Fourier transform of a set of projections along angle θ equals the Fourier transform of the object along the line AB.

FFF

  u and v in frequency domain with w and θ in polar coordinate system, we have: (w, θ)e j2π(x cos θ+y sin θ) w dw dθ (w, θ)e j2π(x cos θ+y sin θ) w dw dθ (w, θ + π)e j2π[x cos(θ+π)+y sin(θ+π)] w dw dθ.

F

  (w, θ + π) = F (-w, θ).

Figure 2 .

 2 Figure 2.15: The bandlimited lter response (a), and the corresponding impulse response (b).

Figure 2 .

 2 Figure 2.16: Fan beam projections where the detectors are positioned with equal space along a straight line.

Figure 2 .

 2 Figure 2.17: Parameters used in the derivation of the reconstruction method for eauispaced fan beam CT imaging.

2 ⃝ 3 ⃝ 2 ⃝; 4 ⃝ 5 ⃝

 23245 Set the reconstruction image an initial value f (0) ; Calculate the forward projection of the f (k-1) (for the rst iteration k = 1, the second iteration k = 2, and so on); Compare the dierence between real measured data and calculated projection by step Correct reconstruction image by back-projection and obtain f (k) ; Increase k by 1 (k = k + 1) and repeat step2

Figure 2 .

 2 Figure 2.18 (b). The second one is dual-source CT technique where two X-ray tubes and two corresponding detectors are equipped in the CT system to simultaneously acquire two set of data for the same slice. These two techniques both use dierent X-ray spectra and their principles are shown in Figure2.18 (b). The third is dual-layer (sandwich) detector technique, which still uses energy-integrating detector, but two sets of scintillation arrays

  Hazirolan et al. applied DECT to the detection of myocardial iron deposition in Thalassaemia patients; the measured myocardial DECT density values showed strong negative correlation with MRI T2* values[START_REF] Hazirolan | Value of dual energy computed tomography for detection of myocardial iron deposition in thalassaemia patients: initial experience[END_REF]. Luo et al. also proved the diagnostic potential of DECT in liver iron content determination using virtual iron content imaging with fty-six patients suspected of having liver iron overload[START_REF] Luo | Dual-energy ct for patients suspected of having liver iron overload: can virtual iron content imaging accurately quantify liver iron content?[END_REF]. However, DECT is unable to give the accurate concentration of iron except for the correlated Hounseld Unit (HU) values and usually requires relatively high dose levels.Dierent from DECT, spectral CT employs photon counting detectors (PCDs) and is able to obtain spectral information of several energy bins with single acquisition, see. Spectral CT technology became available only within the recent years. Thanks to the development of detector technology, PCDs were realized to extract multi-energy information. The rst commercial photon-counting system was introduced for mammography in 2003[Åslund et al., 2007]. Several major manufacturers were also involved in the investigation of energy-resolved CT. In 2013, Philips introduced the IQon Spectral CT, which is the world's rst spectral detector CT. Spectral CT imaging has the potential to overcome many limitations of the previous techniques and apply for many new applications. It can

Figure 2 .

 2 Figure 2.18: Data acquisition from dierent CT systems: conventional CT system where integrating detectors acquire one single sinogram (a); DECT using rapid kilovoltage switching technique or dual-source technique (b) or dual-layer detector technique (c) acquiring two set of projection data; spectral CT system (d) with photon counting detectors acquire several datasets of selected energy bins.

Figure 2 .Figure 2 .

 22 Figure 2.19 illustrates the principle of the PCD. It consists of two parts, a sensor and an electronic part of application-specic integrated circuit (ASIC). The principle of the sensor is the same as the semiconductor based integrating detector that we have discussed in 2.2.3: rstly, the X-rays interact with the semiconductor material, and result in a number of electron-hole pairs; then the electrons and the holes move to the opposite directions, the anode and the cathode respectively; the moving charge carriers induce the signals on the electrodes. The charge signals are then passed to the electronic part (the photon counting ASIC): a preamplier and a shaper convert the charge signals into pulses withamplitudes proportional to the incident photons energies, then the comparators divide the pulses into dierent pre-selected bins according to their amplitudes, the counters record the total number of pulses for each bin and output the counting signal.

2 .

 2 Principles of conventional X-ray CT. The rst CT scanner was invented by Hounseld in 1971. The CT scanner takes measurements of integral projections of an object from dierent angles and then reconstructs the 2-D cross-sectional image through complex mathematics computations. The basic components are X-ray source, object and detector. Currently used energy integrating detectors include scintillator-based and semiconductor-based detector. The most important step in CT imaging is image reconstruction, which retrieves the cross sectional attenuation coecients from the measured X-ray data.

Figure 3 . 1 :

 31 Figure 3.1: Energy dependence of X-ray attenuation coecient for iron (a) and iodine (b).

  above negative log-likelihood function, the line integral of object density A α can be solved. As we can see in Equation(3.5), A α is the line integral of a α , thus a α can be obtained by applying a conventional reconstruction on A α .The K-edge imaging technology enables specic imaging of high-Z K-edge materials, such as iodine bismuth, gold and gadolinium. These materials are good candidates of contrast agents. In the work of Feuerlein, a partially occluded stent was simulated by using a calcied plaque isoattenuated to a surrounding gadolinium chelate solution. By K-edge imaging, the separated gadolinium image clearly shows perfused lumen of vessel[START_REF] Feuerlein | Multienergy photon-counting k-edge imaging: potential for improved luminal depiction in vascular imaging[END_REF]. In the work of Pan et al., brin-specic bismuth-enriched K-edge nanocolloid (nanoK (Bi)) particles are well enhanced in both spectral CT scans of a blood vessel phantom and carotid artery endarterectomy specimens[START_REF] Pan | Computed tomography in color: Nanok-enhanced spectral ct molecular imaging[END_REF].Cormode et al. made in vivo experiments using rabbits after injections of gold nanoparticles and iodinated contrast agent, their results show that gold and iodine distributions are clearly dierentiated[START_REF] Cormode | Multicolor spectral photon-counting computed tomography: in vivo dual contrast imaging with a high count rate scanner[END_REF].

Figure 3 . 2 :

 32 Figure 3.2: Illustration of projection domain decomposition and image domain decomposition. The example here assumes that the spectral CT has 6 energy bins and 6 sinograms are obtained for one acquisition. The former method rstly decomposes acquisition data into M location-dependent target sinograms (M = 4) in this case, and then perform a reconstruction on each obtained sinogram to get the nal spatial image. Image domain method rstly reconstructs the spatial images from their respective projection data corresponding to a given energy bin, and then performs material decomposition on these image domain data (i.e. the reconstructed spatial images).

P

  α (sx, sy) = arg minPα(sx,sy) N ∑ i=1 [ln(λ i (P C α )) -ln(m C i )] 2 .

Figure 3 .Figure 3 . 3 :

 333 Figure 3.3: Illustration of the proj-PR-LLS method for patch C. Patch C is noted by rectangular regions in the acquisition sinograms and the decomposed images. The objective function is calculated by summing up the squared residuals from acquisition data and the regularization term from the target decomposed data. Decomposed images in patch C is obtained by minimizing the objective function.

Figure 3 .Figure 3 . 4 :

 334 Figure 3.4 is an illustration of the ima-PR-LS method. The reconstructed images are divided into location-dependent small patches. At the same time, the decomposed images are also supposed to have the same division. For patch C, the objective function consists of two parts, squared residuals of the reconstructed linear attenuation coecients and regularization term of the decomposed images. Material density within patch C can be obtained by minimizing this objective function. By repeating the decomposition process for all patches, we can get the whole density distribution of all materials.

3. 3 SummaryFigure 3 . 5 :

 335 Figure 3.5: The material decomposition methods discussed and proposed (in blue) in this chapter. The term CR in refers to computed radiography.

Figure 4 .

 4 Figure 4.1 illustrates the basic principle of simulation technique used in VXI. The X-rays are emitted by the source towards the object and pixelized detector. The object consists of dierent small voxels. Each ray interacts with a set of voxels when crossing the object and then the transmitted number of photons is detected by the detector. For the given detector pixel K in Figure4.1, consider the ray SK between source point S and detector K, the expected number of photons N (E) that arrive at pixel K is calculated by Beer-Lambert law:

Figure 4 .

 4 Figure 4.2 illustrates the simulated imaging chain. The users can dene the geometry of the source including shape, size and position, the spectrum and photon ux. Beam parameters enable, if need be, collimators to be taken into account. Object geometry allows import of standard Computer-aided Design (CAD) les and denition of various

Figure 4 . 1 :Figure 4 . 2 :

 4142 Figure 4.1: Illustration of the simulation principles. The ray SK intersects two meshes at points A and B. Geometrical calculations enable the attenuation path length AB to be determined. Ray (1): transmitted photons. Ray (2) and (3): scattered photons. Figureisretrieved from[START_REF] Duvauchelle | A computer code to simulate x-ray imaging techniques[END_REF].

Figure 4 . 3 :

 43 Photon energy (keV)

Figure 4 . 4 :

 44 Figure 4.4: Scheme of sample denition process: geometry of each part of the sample is designed separately, then assembled together to obtain the sample geometry, nally we can dene the material for each part of sample to get the completed sample.

Figure 4 . 5 :

 45 Figure 4.5: Histopathological characteristics of a ruptured plaque (A, B) and the plaque which is prone to disruption (C-D). In images B and D, the plaque area has been colored

Figure 4 .

 4 Figure 4.6 denotes the location of the slice being scanned in our simulation. To build a complete thorax phantom, geometry information is not enough. Material components of dierent parts included in the thorax phantom are of the same importance. The organ and tissue-specic elemental compositions have been published in Publication 44 of International Commission on Radiological Units and Measurements (ICRU44) [ICRU, 1989] and

Figure

  Figure4.7 (a). Table4.1 is a list of this information.

Figure 4 . 6 :

 46 Figure 4.6: Illustration of VHP-Female v2.2 model. Geometries of dierent body parts have been described in this model, including the organs, bones, muscles, tissues, vessels and so on. The blue line in this gure denotes the slice that is scanned during our later simulation.

Figure 4

 4 Figure 4.7: Cross section of the thorax phantom with illustration of dierent material types included (a) and the magnied heart region that contains extra calcium and iodine inserts with dierent concentrations (b).

Figure 4 .

 4 Figure 4.8: Scheme of spectral CT acquisition with computational human phantom.

Figure 4 . 9 :

 49 Figure 4.9: Material decomposition results of proj-LS method (a), proj-LLS method (b)and proj-PR-LLS method (c). Each column represents the same material basis, i.e. soft tissue for the rst column, calcium for the second and iodine for the third.

Figure 4 .

 4 Figure 4.10: Illustration of the ROI and BG regions chosen for quantitative analysis of calcium (a) and iodine (b). The red circles represent ROI and the cyan circles BG regions.

Figure 4 .Figure 4 .

 44 Figure 4.11: Average errors of decomposed densities of calcium (a) and iodine (b) using proj-LS, proj-LLS and proj-PR-LLS methods.

Figure 4 .Figure 4 .

 44 Figure 4.13 illustrates the relationships between calculated average density within each insert and the known true values of calcium and iodine using proj-PR-LLS method. Good linear relationship can be noticed, the slopes of the tting curves are around 1 and the constant intercepts are very small, indicating that the measured densities of calcium and iodine are close to true values.

Figure 4 .

 4 Figure 4.14: Color overlay image of soft tissues(gray), calcium (red) and iodine (green).

Figure 4 .

 4 16 shows the mass attenuation coecients of the components contained in the phantom: ABS, Br, Cl and P, which will be used to distinguish dierent

Figure 4 .Figure 4 .

 44 Figure 4.15: Illustration of the ABS-FR phantom used for spectral CR imaging. Materials of these cubes are introduced inTable 4.2.

Figure 4 .

 4 Figure 4.17: Scheme of spectral CR system geometry. This gure is only a simplied schematic view of the simulation system where the dimensional information (D SC , D DC , size of phantom and detector) is not in realistic scale.

Figure 4 .

 4 18 (a). It is observed that all cubes are visible in the ABS basis image, due to the presence of ABS in either mixtures (ABS-FR) or pure ABS. The ABS-TBBPA cubes are well separated and highlighted in the Br basis image. When cube thickness d = 0.3 and 0.5 mm, the ABS-TBBPA cubes and ABS-DDC-CO tubes are visible in the Cl basis image; When cube thickness increases to 1 mm and 2 mm, all three ABS-FR cubes are visible, however, the signal intensity of ABS-TBBPA cubes in the Cl basis image decreases compared to those with smaller thickness; When d = 4 mm, only ABS-DDC-CO cubes and ABS-RDP cubes are visible in the Cl basis image. To conclude, in noise-free situation, ABS and Br can be well identied; the determination performance of other FRs improves with object thickness: when d = 4 mm, the FRs containing Cl and P are present in the Cl basis image without cross-talk of ABS-TBBPA material, one can further distinguish the two FRs according to their large density dierence.

Figure 4 .

 4 Figure 4.18: Decomposition results under dierent conditions using proj-PR-LLS method. (a): simulation of CR without noise, r = 0; (b): simulation of CR with Poisson noise, r = 0; (c): simulation of CR with Poisson noise, r = 10 4.4 . The rst column represents ABS basis, second column Br basis, and third column Cl basis image. The dash lines on images of (a) mark the position where 1-D proles are plotted in later analysis.

Figure 4 .

 4 Figure 4.19: Flowchart of the decomposition method with classication for one patch. fval1 and fval2 represent the cost function values obtained at the end of decomposition 1 and decomposition 2.

Figure 4 .

 4 Figure 4.19 shows the procedure of this method to obtain values for one patch. Instead of the 3-material (ABS, Br and Cl) decomposition using proj-PR-LLS method, the rst

Figure 4 .

 4 Figure 4.20 shows the decomposition results of proj-PR-LLS-clas method with the same condition as proj-PR-LLS method, whose results has been shown in Figure4.18.In noise-free condition (Figure4.20 (a)), ABS and Br contained FR can be identied from the ABS and Br basis images. Meanwhile, the cubes of FR materials containing Cl and P appear in the Cl basis image, but they have signicant density dierence for cubes with the same thickness, therefore, they can be easily distinguished by the observer even if the concentration of FRs or thickness of cubes changes more or less. Figure4.21 gives the 1-D proles along dash lines in Figure4.18 (a) and4.20 (a), the comparison between theoretical density integrals (P α in Table4.3) of 3 basis materials and our calculations

Figure 4 .

 4 20 (b) that the ABS-TBBPA cubes are indicated in the Br basis image, with some mistakenly appeared pixels of the ABS-DDC-CO cubes. This situation is improved obviously when we use r = 10 4.4 : less pixels of ABS-DDC-CO cubes appears in the Br basis image of Figure 4.20 (c). Besides, comparing the Cl basis images of (b) and (c), (c) is closer to the Cl basis image of (a), especially when cube thickness d = 2 mm and 4 mm where only the ABS-DDC-CO cubes and ABS-RDP cubes are highlighted. Therefore, a good choice of relaxation parameter can improve the decomposition performance.

Figure 4 .

 4 Figure 4.20: Decomposition results under dierent conditions using proj-PR-LLS-clas method. (a): simulation of CR without noise, r = 0; (b): simulation of CR with Poisson noise, r = 0; (c): simulation of CR with Poisson noise, r = 10 4.4 . The rst column represents ABS basis, second column Br basis, and third column Cl basis image. The dash lines on images of (a) mark the position where 1-D proles are plotted in later analysis.

Figure 4 .

 4 Figure 4.21: Performance comparison of proj-PR-LLS method and proj-PR-LLS-clas method in noise-free condition: 1-D proles along the dash lines in Figure 4.18 and 4.20.ABS basis (top), Br basis (middle) and Cl basis (bottom). Black curves represent the theoretical density integrals of basis materials (P α for thickness of 1 mm in Table4.3) if

Figure 5 . 1 :

 51 Figure 5.1: Illustration of spectral CT acquisition with 4 (a) and 6 (b) energy bins for the same X-ray spectrum.

Figure 5 . 2 :

 52 Figure 5.2: Phantom illustration. Each row contains the same material or materials soluted by water at dierent concentrations, types of material are indicated on the left hand and concentrations are annotated inside each insert on the gure with unit of mg/cc. Symbol `#' marks the mixture inserts which combine iron, calcium, potassium and water together.

Figure 5 .Figure 5 . 3 :

 553 Figure 5.3 shows the mass attenuation coecients of dierent materials contained in this phantom within the detection energy range, i.e. 30 keV-90keV. It can be observed that all materials are without K-edge, and those with close atom numbers like Ca and K, water and PMMA have similar attenuation properties, making them dicult to be separated from each other. In this case, the attenuation caused by either of the component is supposed to be represented by the same basis material.

•Figure 5 . 4 :Figure 5 . 5 :

 5455 Figure 5.4: Scheme of spectral CT acquisition.

Figure 5 . 7 :

 57 Figure 5.7: Comparison of the decomposed average density with true values by proj-PR-LLS method: iron (a) and potassium (b) when acquisition parameters are 6 energy bins, 1 mA•s per projection; iron (c) and potassium (d) when acquisition parameters are 60 energy bins, 10 mA•s per projection. The error bars indicate the standard deviation of decomposed density.

Figure 5 .Figure 5

 55 Figure 5.8 is a comprehensive illustration of the average errors of identied iron (in the mixture) concentrations with all dierent acquisition parameters. Quantication error tends to decrease for iron of larger concentration. Optimal performance can be observed for 200 mg/cc iron, with exposure factor of 100 mA•s, which means that the improvement of quantication thanks to a higher exposure factor has a major inuence, even if the number of energy bins increases.

Figure 5 . 9 :

 59 Figure 5.9: CNRs of decomposed iron images for 100 mg/cc of iron with dierent acquisitions using proj-PR-LLS method. Each group of acquisitions has the same number of energy bins but dierent exposure factors.

Figure 5 .

 5 Figure 5.9 plots the CNRs of the decomposed iron images for 100 mg/cc of iron. Clearly

Figure 5 .Figure 5 .Figure 5 .

 555 Figure 5.10: Decomposed density in unit of mg/cc of water basis (left column), iron basis (middle column) and potassium basis (right column) by ima-PR-LS method. Data acquisition parameters: 6 energy bins, 1 mA•s per projection (a), and 60 energy bins, 10 mA•s per projection (b). The red, blue and magenta lines on the iron basis image mark the position where 1-D proles are plotted in later analysis.

Figure 5 .Figure 5 .Figure 5 .

 555 Figure 5.10 (b) shows the decomposition results when N =60, and F =10 mA•s. Obvious dierence can be observed when comparing withFigure 5.10 (a). Only two rows of inserts exist in the iron basis image (middle) of Figure5.10 (b), representing iron in the mixture and in water; for the potassium basis image (right), inserts containing calcium and potassium are highlighted. Figure5.13 gives the relationship between ρins of Fe/K and the true values. Similarly to those obtained by proj-PR-LLS method in Figure5.7, good linear properties can be observed; for iron soluted in the mixture or in water, the slopes of tting curves were almost 1, and the intercepts were considerably small; for potassium, the calculated concentrations in mixture inserts were almost two times the values in water inserts, caused by the fact that calcium was also classied into potassium basis.

Figure 5 .Figure 5 .

 55 Figure 5.16 demonstrates the computational time of proj-PR-LLS and ima-PR-LS methods for dierent acquisitions. The computation time is measured on a computer with processor Intel Core i5 and 8 GB RAM. The software platform used was Matlab R2017a. It can be observed that the computation time of proj-PR-LLS method depends largely on the number of energy bins while the exposure factor has little inuence. With the same acquisition parameters (N = 60 and the same F ,), ima-PR-LS method takes

Figure 5 .

 5 Figure 5.17: Zoomed views of the decomposed iron images (N = 60, F =10 mA•s) obtained by proj-PR-LLS (a) and ima-PR-LS (b) methods to show the pattern defects caused by patchwise regularization.

  (a)) since the image is obtained by applying FBP reconstruction to the decomposed sinogram.

6. 1

 1 ConclusionIn this thesis, we have studied the quantitative material decomposition methods for spectral CT. As mentioned in Chapter 1, the main contributions of the present work lie in three aspects: (1) study of the material decomposition methods for spectral X-ray CT;(2) validation and assessment of the proposed methods via simulation for medical and industrial applications; (3) investigation of the inuence of dierent acquisition parameters, i.e. exposure factor and number of energy bins, on material decomposition performance.Before presenting the achievements of the present thesis work, the basic principles of spectral CT have been reviewed in Chapter 2. Based on the fact that materials have dierent attenuation properties to X-ray, CT scanners take measurements of integral projections of an object from dierent angles and then reconstruct the 2-D cross-sectional image through more or less complex computations. Conventional CT devices are equipped with energy integrating detectors which provide integrated information of all the transmitted beam while spectral CT employs photon counting detectors and is able to obtain spectral information of several energy bins with single acquisition. Spectral CT enables us to overcome many limitations of the conventional CT techniques and opens up many new application possibilities, among which quantitative material decomposition is the hottest topic.

FigureFigure A. 1 :

 1 Figure A.1: Illustration of ROI for the determination of r.

α

  Figure A.2: The L-curve with r varying from 10 -10 to 10 10 (a), and from 10 4 to 10 5 (b).

  Différents paramètres d'acquisition, c'est-à-dire le facteur d'exposition et le nombre de bandes d'énergie, ont été simulés pour étudier leur influence sur la performance de décomposition pour la détermination du fer. MOTS-CLÉS : Tomographie, Rayons X, CT spectrale, Détecteur à comptage des photons, Décomposition de matériaux, Simulation Laboratoire (s) de recherche : Laboratoire Vibrations et Acoustique (LVA) Directeur de thèse: KAFTANDJIAN Valérie Président de jury : DESVIGNES Michel Composition du jury : LAQUERRIERE Patrice, VINCENT Nicole, DOUEK Philippe, DESVIGNES Michel, RUAN Su, KAFTANDJIAN Valérie, DUVAUCHELLE Philippe, ZHU Yuemin

  extrémités et de la tête est de 350 à 450[START_REF] Zhao | Noise, sampling, and the number of projections in cone-beam ct with a at-panel detector[END_REF]. Dans la présente étude, 1200 projections ont été utilisées dans le balayage du fantôme simulé pour une meilleure qualité d'image. Ce nombre pourrait conduire à une dose de patient interdite et ceci devrait être considéré dans de futures études.Nous avons simulé un fantôme composé de plusieurs cubes d'une hauteur de 10 mm, d'une largeur de 10 mm et d'épaisseurs diérentes. Comme le montre la Figure10, chaque colonne de cubes est du même matériau (noté en haut) et chaque rangée de cubes est de même épaisseur (notée sur la gauche).

	3.2 Application industrielle: tri plastique ABS-RF ABS-ABS-ABS-
	Les déchets d'équipements électriques et électroniques (DEEE) ont augmenté rapidement TBBPA (Br) DDC-CO (Cl) RDP (P) ABS
	en raison du développement de l'industrie électronique. Ces déchets ont diverses compo-
	sitions parmi lesquelles les polymères jouent un rôle important. L'acrylonitrile-butadiène-Thickness
	d=0.3mm
	styrène (ABS) est un composant majeur parmi tous les polymères. En raison de la présence
	d=0.5mm
	d=1mm
	d=2mm
	d=4mm	10 mm
	10 mm

Trois images ont été obtenues en utilisant l'une ou l'autre méthode: image des tissus mous, du calcium et de l'iode. Les résultats ont montré que toutes les méthodes peuvent donner une distribution quantitative des trois matériaux, parmi lesquels l'image des tissus mous contient de la peau, du tissu adipeux, du tissu mou moyen, du cartilage, de l'aorte, du poumon et du coeur; l'image du calcium met bien en évidence les os et les solutions de calcium insérées; L'image de l'iode contient seulement les solutions d'iode insérées. L'analyse quantitative, incluant les erreurs moyennes et les CNR, a démontré que la méthode proj-LLS proposée et la méthode proj-PR-LLS ont une meilleure précision quantitative et une meilleure qualité d'image que la méthode proj-LS couramment utilisée. De plus, la méthode proj-PR-LLS a la capacité d'améliorer les CNR des images décomposées tout en obtenant une meilleure précision quantitative (calcium) ou inférieure mais comparative (iode), en comparaison avec la méthode proj-LLS. Il y a des améliorations potentielles à notre travail. La dose de rayonnement pour le patient est liée au nombre de projections. Un nombre inférieur de projections peut réduire la dose pour le patient mais également diminuer la qualité de l'image en raison de l'échantillonnage plus clairsemé des données d'imagerie [Van Daatselaar et al., 2004]. Le nombre typique de projections utilisées en CT à faisceau conique dans les applications du sein, des d'énergie des EEE, les matériaux en plastique doivent répondre à des normes de sécurité incendie élevées. Cependant, il n'est pas possible de les réaliser avec des polymères purs, donc des retardateurs de ammes (RF) sont ajoutés pour modier l'inammabilité des plastiques et augmenter la résistance au feu. Le recyclage des plastiques issus des DEEE est dicile en raison de l'existence des RF (en particulier les RF bromés et chlorés), ce qui pourrait entraîner une pollution environnementale grave. Par conséquent, la technique de tri est importante pour le traitement ultérieur et le recyclage des DEEE an d'éviter les problèmes environnementaux et de santé majeurs.

Dans le présent travail, nous allons étudier la capacité de décomposition en matériaux de la méthode proj-PR-LLS pour l'imagerie CR spectrale dans l'application de tri des plastiques ABS-RF (c'est-à-dire en radiographie et non pas en tomographie puisque les déchets délent sur un convoyeur à grande vitesse). En outre, nous développons également une méthode de décomposition basée sur la classication pour cette application spécique, appelée méthode proj-PR-LLS-clas. Ses performances seront évaluées par rapport à la méthode proj-PR-LLS.

3.2.1 Fantôme simulé & géométrie du système

Comme nous l'avons dit dans la section ci-dessus, l'ABS est un composant majeur parmi tous les polymères utilisés dans les équipements électriques et électroniques, nous l'avons donc choisi comme matériau polymère pour l'investigation. Nous avons également choisi trois types de retardateurs de amme couramment utilisés, y compris les ignifugeants bromés (Br), chlorés (Cl) et phosphorés (P). Ce sont respectivement le tétrabromobisphénol A (TBBPA), le déchlorane plus (DDC-CO) et le bis diphénylphosphate de résorcinol (RDP). Trois matériaux ABS-RF simulés ont été obtenus en mélangeant chaque RF avec de l'ABS à un pourcentage massique de 15%. Ces choix sont issus d'une étude précédente de CT en bi-énergie menée au laboratoire lors de la thèse de Florian Montagner.

6 bins 10 bins 15 bins 20 bins 30 bins 60 bins Different acquisitions 10 -1 10 0 10 1 10 2 10 3 CNR

	Figure 19: CNR des images de fer décomposées pour 100 mg/cc de fer avec diérentes
	acquisitions utilisant la méthode proj-PR-LLS. Chaque groupe d'acquisitions a le même
	nombre de bandes d'énergie mais des facteurs d'exposition diérents.
	4.2.2 Méthode ima-PR-LS
	Pour la méthode ima-PR-LS, nous choisissons l'eau, le potassium et le fer comme matériaux
	de base. Le paramètre de relaxation r est déni à 1 et la taille de chaque zone pour la
	régularisation est 5 × 5. Avec divers paramètres d'acquisition, les résultats de la simulation
	montrent qu'avec la méthode de décomposition ima-PR-LS, le fer ne peut être séparé du
	calcium et du potassium que lorsque le nombre de bandes est 60. Nous montrons sur la

  dx dy, δ is the Dirac delta function. P θ (t) is called the Radon transform of function f (x, y). P θ (t) is physically the projection of f (x, y)

	(2.17)

where

  along the line which is in direction θ and has distance t from the origin of coordinate. If we x θ and change t to all dierent values, we will obtain a set of projection data in direction θ; If we change θ at the same time we will obtain sets of projection data from dierent directions, that is called sinogram.

	Projection data
	P θ (t)	y
		θ
		x
	µ(x, y) or
	f(x, y)

Table 2 .

 2 1 gives a comparison of several PCDs.

Table 2 .

 2 1: Comparison of dierent PCDs. This table is not a complete list of all available PCDs and is only for comprehensive review.

	Energy	resolution	NA	4.75% (5.8 keV	FWHM @	122 KeV)	7 keV	FWHM @	60 KeV	NA	12%-26%	NA	≤15%	FWHM @	122 KeV
	Maximum	count rate	(Mcps/pixel)	1-2		6			20		NA	NA	0.145-0.21		2
	No. of	energy	thresholds	6		2			4		4	2	1, 2, 4, 8		5
	Index Name Manufacturer Material Pixel number Pixel size (µm×µm)	1 MEXC (Refs group 1) 1 Gamma Medica Ideas, Northridge, California CdTe 1×1024 380×160	2 DXMCT-1 (Refs group 2) 2 DxRay Inc, Northridge, California CdTe 16×16 1000×1000	3 DXMCT-2 (Refs group 3) 3 DxRay Inc, Northridge, California CZT 16×16 500×500	4 Siemens 2010 (Refs group 4) 4 Siemens AG, Forchheim, Germany CdTe 64×64 225×225	5 MicroDose SI (Refs group 5) 5 Sectra Mamea AB, Solna, Sweden Si strip NA 50×50	6 Medipix3 (Refs group 6) 6 CERN, Geneve, Switzerland CMOS 256×256 55×55 110×110	7 C10413 (Refs group 7) 7 Hamamatsu Photonics, Japan CdTe 1×64 800×500	1 Reference 2 group 1 [Feuerlein et al., 2008, Schlomka et al., 2008, Pan et al., 2012] Reference group 2 [Barber et al., 2009, Iwanczyk et al., 2009] 3 Reference group 3 [Barber et al., 2013]

4

Reference group 4

[START_REF] Kappler | A research prototype system for quantum-counting clinical ct[END_REF], Kappler et al., 2012]

5

  Iterative reconstruction has the potential to reduce noise and improve image quality, but it is not applicable for widely used practical applications due to huge computation costs.

	Making use of multi-energy information, spectral CT has the potential to discriminate
	3. Principles of spectral X-ray CT. Spectral CT employs photon counting detectors and dierent components inside the object. Many material decomposition methods have been
	is able to obtain spectral information of several energy bins with single acquisition. proposed and dierent experimental systems are under development for spectral CT. Ac-
	This technology became available only within the recent years thanks to the develop-ment of detector technology. PCDs are semiconductor based detectors that have the +D=FJAH ! cording to the formulation of the decomposition, there are mainly two kinds of methods.
	ability to discriminate transmitted photons corresponding to selected energy bins.
	Spectral CT imaging enables to overcome many limitations of the previous tech-niques and opens up many new applications. Quantitative material decomposition is the most investigated topic for spectral CT, including K-edge imaging of con-Material decomposition methods
	trast materials and selective imaging of other previously indistinguishable materials.
	Other applications include energy weighting and virtual monochromatic imaging to
	improve image quality.

There are two main categories of methods, analytical reconstruction (FBP reconstruction) and iterative reconstruction. Page 35 Cette thèse est accessible à l'adresse : http://theses.insa-lyon.fr/publication/2018LYSEI056/these.pdf © [T. Su], [2018], INSA Lyon, tous droits réservés Contents 3.1 State of the art . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.1.1 Eect-based decomposition / K-edge imaging . . . . . . . . . . 3.1.2 Material-based decomposition . . . . . . . . . . . . . . . . . . . . 3.2 Patchwise regularized decomposition . . . . . . . . . . . . . . . .

3.2.1 Projection domain . . . . . . . . . . . . . . . . . . . . . . . . . . 3.2.2 Image domain . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.3 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

Table 4 .

 4

	Body tissue	H	Elemental composition (% by mass) C N O Others	Density (kg/m 3 )
	Lung	10.3	10.5	3.1	74.9	Na(0.2); P(0.2); S(0.3); Cl(0.3); K(0.2)	1050
	Heart	10.3	12.1	3.2	73.4	Na(0.1); P(0.1); S(0.2); Cl(0.3); K(0.2); Fe(0.1)	1060
	Average soft tissue	10.6	30.8	2.4	55.5	Na(0.1); P(0.1); S(0.2); Cl(0.1); K(0.2)	1010
	Adipose	11.4	59.8	0.7	27.8	Na(0.1); S(0.1); Cl(0.1)	950
	Skin	10	20.4	4.2	64.5	Na (0.2); P(0.1); S(0.2); Cl(0.3); K(0.1)	1090
						Mg(0.1);P(4.0);S(0.2);
	Sternum	7.8	31.6	3.7	43.8	Cl(0.1); K(0.1); Ca(8.5);	1250
						Fe(0.1)
	Cartilage	9.6	9.9	2.2	74.4	Na(0.5); P(2.2); S(0.9); Cl(0.3)	1100
						Na(0.1);Mg(0.1); P(6.0);
	Ribs	6.4	26.3	3.9	43.6	S(0.3); Cl(0.1); K(0.1);	1410
						Ca(13.1)
	Aorta (blood)	10.2	11	3.3	74.5	Na(0.1); P(0.1); S(0.2); Cl(0.3); K(0.2); Fe(0.1)	1060
						Na(0.1); Mg(0.1); P(6.1);
	Vertebra	6.3	26.1	3.9	43.6	S(0.

1: Elemental compositions and mass densities of body tissues included in the thorax phantom.

Table 4 .

 4 2: Components of three ABS-FR materials used for the phantom.

		Material(ABS-FR)	ABS-TBBPA	ABS-DDC-CO	ABS-RDP
	Density of material ρ(mixture)	1060 mg/cm 3	1060 mg/cm 3	1060 mg/cm 3
		Chemical formula of FR	C 15 H 12 Br 4 O 2	C 18 H 12 Cl 12	C 30 H 24 O 8 P 2
		Mass % (ω)of FR	15%	15%	15%
	Mass % (ω) of Br, Cl and P, respectively	8.82%	9.76%	1.62%
	ρ eff	* of Br, Cl and P, respectively	93.5 mg/cm 3	103.5 mg/cm 3	17.2 mg/cm 3

* : ρ eff is the eective density calculated by ρ eff (α) = ρ(mixture) × ω(α).

Table 4 .

 4 3: Theoretical values of P α (unit: mg/cm 2 ) for each material. Material ABS has two columns of values, the rst is ABS pure, which corresponds to the material of the 4th column of cubes in the phantom shown in Figure 4.15, and ABS in mixture corresponds to the ABS components in the other cubes of Figure 4.15.

	Thickness	ABS (Pure)	ABS in mixtures	Br in ABS-TBBPA	Cl in ABS-DDC-CO	P in ABS-RDP
	0.3 mm	31.80	27.03	2.81	3.11	0.52
	0.5 mm	53.00	45.05	4.68	5.18	0.86
	1 mm	106.00	90.10	9.35	10.35	1.72
	2 mm	212.00	180.20	18.70	20.70	3.44
	4 mm	424.00	360.40	37.40	41.40	6.88

  . With the objective of quantifying iron from other materials, we simulate the acquisition process with dierent parameters, i.e. exposure factor and number of energy bins, and investigate their inuence on the decomposition performance.Iron is a micronutrient essential for adequate erythropoietic function, oxidative metabolism and cellular immune response to bacterial infection. Disorders in iron, such as iron deciency or overload, are important risk factors for heart diseases such as coronary artery disease and heart failure[Von Haehling et al., 2015, Muñoz et al., 2011]. There are dierent methods for the evaluation of tissue iron content, which can be classied into two cate-

	5.2 Method and material
	5.2.1 Medical context

Page 84 Cette thèse est accessible à l'adresse : http://theses.insa-lyon.fr/publication/2018LYSEI056/these.pdf © [T. Su], [2018], INSA Lyon, tous droits réservés solutionsgories: invasive techniques and noninvasive techniques. Invasive techniques include biopsy and serum ferritin level. Noninvasive techniques include superconducting quantum interference device (SQUID) and MRI

  Hazirolan et al. applied DECT to the detection of myocardial iron deposition in Thalassaemia patients; the measured myocardial DECT density values showed strong negative correlation with MRI T2* values[START_REF] Hazirolan | Value of dual energy computed tomography for detection of myocardial iron deposition in thalassaemia patients: initial experience[END_REF]. Luo et al. also proved the diagnostic potential of DECT in LIC determination using virtual iron content imaging with fty-six patients suspected of having liver iron overload[START_REF] Luo | Dual-energy ct for patients suspected of having liver iron overload: can virtual iron content imaging accurately quantify liver iron content?[END_REF]. However, DECT is unable to give the accurate concentration of iron except for the correlated Hounseld Unit (HU) values and usually requires relatively high dose levels. .2. METHOD AND MATERIAL the carotid atherosclerotic plaque; water-like, calcium-like and lipid-like tissues were well separated according to histology comparison, but iron was not distinguished from calcium[Zainon, 2012].Li et al. developed an image domain material decomposition algorithm

	Unlike DECT, spectral CT employs PCD and is able to obtain spectral information
	of several energy bins with single acquisition. Its ability in material decomposition has
	been discussed in previous chapters, including K-edge material and non-K-edge material
	determination. Iron, as a non-K-edge material within diagnostic X-ray energy range (10

keV-150keV), is more challenging to be distinguished.

Zainon et al. used

small animal spectral micro-CT with Medipix3 detector (thresholds: 10, 16, 22 and 28keV) to scan Page 85 Cette thèse est accessible à l'adresse : http://theses.insa-lyon.fr/publication/2018LYSEI056/these.pdf © [T. Su], [2018], INSA Lyon, tous droits réservés

5

•

  Ting Su, Valérie Kaftandjian, Philippe Duvauchelle, and Yuemin Zhu. A spectral Xray CT simulation study for quantitative determination of iron. Nuclear Instruments and Methods in Physics Research Section: A, 2018, 894:39-46 (Published) • Bingqing Xie, Ting Su, Valérie Kaftandjian, Pei Niu, Feng Yang, Marc Robini, Yuemin Zhu, Philippe Duvauchelle. Material Decomposition in X-ray Spectral CT Using Multiple Constraints in Image Domain. Journal of Nondestructive Evalua-Ting Su, Valérie Kaftandjian, Philippe Duvauchelle, Philippe Douek, and Yuemin Zhu. Quantitative material decomposition method for spectral CT imaging. 7th Conference on Industrial Computed Tomography (iCT 2017), Leuven, Belgium, Feb 2017. (Published) • Ting Su, Valérie Kaftandjian, Philippe Duvauchelle, and Yuemin Zhu. Material decomposition for spectral X-ray CT: application to calcium and iodine identication. 2017 IEEE Nuclear Science Symposium and Medical Imaging Conference (IEEE NSS/NIC 2017), Atlanta, USA, Oct 2017. (Published)• Ting Su, Valérie Kaftandjian, Philippe Duvauchelle, and Yuemin Zhu. Plastic sort-

	tion. (Under review)
	Conference papers:

•

ing by X-ray radioscopy with photon counting detector. 12th European Conference on Non-Destructive Testing (ECNDT 2018),

Gothenburg, Sweden, June, 2018 (Accepted)

Cette thèse est accessible à l'adresse : http://theses.insa-lyon.fr/publication/2018LYSEI056/these.pdf © [T. Su],[2018], INSA Lyon, tous droits réservés

Acknowledgements

To make a quantitative analysis, we calculated the average density ρins of each mea- sured iron and potassium insert, the relationship between ρins and true values is shown in Figure 5.7. The corresponding error bar for each average measurement indicates the variability of data. The tting curves of ρins at dierent concentrations against the the- oretical values t ins showed very good linear properties. For iron (Figure 5.7 (a) and (c)), no matter when it was soluted in the mixture or in water, the slopes of tting curves were almost 1, and the intercepts were considerably small, meaning that the quantication was quite accurate. For potassium ( The red diamonds in Figure 5.14 represent the average errors of the identied iron concentration with acquisitions of 60 energy bins and various exposure factors using ima-PR-LS method. Quantication accuracy was improved with the increase of the factor of exposure and iron density. With the acquisition of 100mA•s per projection, the average error of the estimated iron content varies from 35.6% to 9.0% for the concentrations from 15 mg/cc to 200 mg/cc. The CNR of the decomposed iron images for 100 mg/cc of iron varies from 14.1 to 326.0 when increasing exposure in the order, while the CNRs obtained for the same acquisition using proj-PR-LLS method vary from 1.0 to 35.0, as shown in Figure 5.15.

We can observe that with 60 energy bins, ima-PR-LS method exhibits comparable average errors with proj-PR-LLS method, but with much higher CNRs.