Multielement imaging in computerised X-ray tomography

Multielement imaging in computerised X-ray tomography

Nuclear Instruments and Methods in Physics Research A271 (1988) 671-677 North-Holland, Amsterdam 671 MULTIELEMENT IMAGING IN COMPUTERISED X-RAY TOMO...

517KB Sizes 0 Downloads 71 Views

Nuclear Instruments and Methods in Physics Research A271 (1988) 671-677 North-Holland, Amsterdam

671

MULTIELEMENT IMAGING IN COMPUTERISED X-RAY TOMOGRAPHY Joseph FRYAR, Kieran J. McCARTHY and Adrian FENELON School of Physical Sciences, National Institute for Higher Education, Dublin, Republic of Ireland

Received 11 November 1987 A method by which several elements can be independently imaged in a single X-ray CT scan is described. The method is based on measurements of differential absorption across absorption edges. The images produced show the distribution of element concentration. Mass attenuation coefficients for the elements comprising the matrix are not required because corrections for matrix effects are made. The method is capable of imaging elements which are consecutive in the periodic table. Detection limits, in the apparatus described, are of the order of mg/cm3 . The elements chosen for imaging were palladium, silver and cadmium. 1. Introduction The ability of X-ray computerised tomography to produce images which show the spatial variation of linear attenuation coefficient within a specimen is well established . Indeed it is the main use of the technique in medical and industrial imaging. It can, however, give more detailed information about the specimen if measurements are taken at more than one X-ray beam energy or X-ray tube voltage. By changing the tube voltage it is possible to image separately low and high atomic number elements [1,2] but it is usually not possible to distinguish between elements which are close in atomic number . A limited amount of elemental discrimination may be obtained by switching the mean energy of the beam across the K absorption edge of an analyte [3] but this works only if there are no other elements with K edge between the mean beam energies . Glantschnig and Holliday [4] describe a method called "mass fraction profiling" which they applied to porous silica soot boules containing germanium. This method depends upon an accurate knowledge of all mass attenuation coefficients associated with the component parts of the object and does not work well for elements which are close in atomic number. Good elemental discrimination can be obtained if, instead of broad spectrum beams, monochromatic beams are used [5-8]. The beam energies straddle the K (or L) edge energy of the analyte. In a previous paper [7] we described an experiment in which two elements (Pd and Cs) within the same specimen were independently imaged. Nearly monochromatic radiations from radioactive sources were used . Because it was not possible to get the energies of the beams very close to the K edge of the analyte, corrections for matrix effects were applied to the data. In these methods knowledge of the absorp0168-9002/88/$03 .50 © Elsevier Science Publishers B.V . (North-Holland Physics Publishing Division)

tion coefficients of the matrix elements are not required . One of the methods discussed involved extrapolation of data taken close to the absorption edge to the absorption edge energy . It was suggested that a continuous spectrum could be used and it was mentioned that such an experiment would be carried out. The present paper describes the results of this experiment in which three consecutive elements in the periodic table were imaged in a single scan . The resultant tomographs show the spatial variation of the concentration of the analyte elements . Sensitivities of the order of mg/cm3 are possible.

2. Theory If a specimen consists of an analyte element in a matrix material it is possible to produce a CT image, which shows the location of the analyte element, by scanning the object with two monochromatic X-ray beams [5-7]. The energies of the two beams straddle an absorption edge of the analyte. We showed in ref. [7] that the image shows the concentration of the analyte provided the projection data is the equivalent thickness to of the analyte, given by eq . (7) of ref. [7] as ta

_

~ In (

Ni Noh

\ N h N01 ~~

uah - ual)

t

tm( umh - uml) uah-ual

where N is the detected photon flux, a is the mass attenuation coefficient, tm is the equivalent thickness of the matrix and subscripts h and 1 refer to energies above and below the absorption edge respectively. Subscript 0 denotes an incident photon flux before passing through the specimen and m and a denote values appropriate to

67 2

J. Fryar et aL / Multielement imaging

the matrix and analyte respectively . Eq . (1) applies to a pencil beam of X-rays . Eq . (1) is the equation by which the equivalent thickness of the analyte should be calculated from the experimental data . It is not practicable, except possibly in some very limited circumstances, to use this equation . If the mass attenuation coefficients of the analyte are known the first term in the equation is easy to evaluate . However, the second term depends on the equivalent thickness of the matrix and its mass attenuation coefficient. These will in most cases be unknown. Corrections will have to be made to the basic scan data to allow for these matrix effects. Two methods by which matrix effects can be eliminated were discussed in the previous paper [7]. One method, which was actually used in those measurements, is to scan the matrix alone and then to scan with the analyte present. This requires two scans of the object . The equivalent thickness, corrected for matrix effects, is given by eq . (15) of ref. [71 as

Extrapolation to absorption edges. Step 21 Rotation 4

In(No/N)

X-ray energy (keV)

Fig. 1. Variation of ln(NO/N) as a function of X-ray energy near the K absorption edges of palladium, silver and cadmium. N,, is the incident number and N the transmitted number of photons in a region of interest . This is data taken in the regions of interest 1-8 inclusive in fig. 4 for the 21st beam on the 4th rotation . The extrapolated values are used in eq . (7) to de termine the equivalent thickness of the analytes in this particular beam direction.

(2) - In{( NI/Nh)(Nmh/Nmi)) /( 14 ah - U..),

in which Nmh and N.,I are transmitted fluxes through the matrix alone. Such a method is possible only if the analyte can be added to the matrix after the first scan . The second method can be implemented in a single scan . We shall follow the argument in ref. [7]. If a monochromatic X-ray beam is passed through the object the transmitted flux is given by N=Noexp { -(tau .+tmum) ,

(3)

ln(NO/N) = tau . + tmum .

(4)

or

The variation of this factor as a function of energy is shown in fig. 1 in the vicinity of the K absorption edges of palladium, silver, and cadmium. The nature of the matrix used in these measurements will be described later. The results have been extrapolated to the K absorption edges to yield, for each element, the extrapolated values, YI

=

t.ual

(5)

+ tniunal,

and

Yh = ' a Uah + 'muni), .

(6)

In these equations the Y's are extrapolated values of the logarithms and the mass attenuation coefficients are the values extrapolated to the high and low side of the K edge . Across the K edge umh = uml . Thus subtracting eqs. (5) and (6) gives ta_(Yh-

YI)/( uah

-t al) -

(7)

Thus projection data, corrected for matrix effects, can be obtained by scanning the object with a continuous spectrum of X-rays and measuring, at each X-ray beam position, the transmitted photon flux in two narrow

energy bands on either side of the absorption edge of the analyte . ln(NO /N), taken in each energy band, is extrapolated to the absorption edge and then eq . (7) is used to determine the equivalent thickness of the analyte which is the projection data as shown below. In the above it has been assumed that the photon counting apparatus does not suffer from dead-time effects. At high count rates, where it might be expected that dead-times would be important, there is no reason to suppose that the dead time in N is the same as the dead-time in No . It is now shown that dead-time corrections are not required . Suppose that the corrected value of N and No are cN and CNo respectively . The correction factor will be the same for each point on the spectrum . The logarithm in eq . (4) becomes ln(CN0/cN) or ln(C/c) + ln(NO /N). Thus the effect of the dead-time correction is an equal displacement along the ordinate axis for each point in fig. 1 : the extrapolated value of Yh - YI is unaffected . However, high dead-times will result in poorer counting statistics . At any X-ray beam position and analyte equivalent thickness is given by ta=

C(r) dy,

where C(r) is the analyte concentration at the position r and integration is along the line joining the source to the detector . Combining eq . (7) with eq. (8) gives (

Yh - YI )/( Uah

- Ual ) =

C(r) dy .

9

By using the equivalent thickness as the projection data the reconstructed image shows analyte concentration .

67 3

J. Fryar et al. / Multielement imaging

3. Apparatus 3.1 . The specimen The aim of the present series of measurements was to image separately three consecutive elements in the periodic table. The elements chosen were ,Pd, 4,Ag and ,Cd each at two different concentrations . They were mixed with epoxy glue and individually inserted into six test tubes each 5 mm in diameter . These tubes were then inserted parallel to each other into a 30 mm diameter cork. The object of the experiment was to produce cross-sectional images perpendicular to the axes of the test tubes. A separate image was to be produced for each element with 1 mm spatial resolution . 3.2 . The X-ray source and beam filter The K absorption edges of the three elements were 24 .35 keV for palladium, 25 .52 keV for silver and 26 .71 keV for cadmium. The range in energy of the X-rays required was ideally from about 23 to 28 keV, i.e. from about 1000 eV below the K edge of palladium to about 1000 eV above the K edge of cadmium. The X-ray source used was a Philips type PW2184/00 tube . The maximum operating voltage on this tube is 60 kV with tube currents up to a few tens of milliamperes . The target was tungsten and the exit window was 1 mm thick beryllium. Any X-ray photons with energies outside the range 23-28 keV were superfluous to the experimental requirements . Thus the X-ray beam was heavily filtered. Tin has a K absorption edge energy of 29 .19 keV. The variation of the mass attenuation coefficient of tin in the vicinity of the K edge is shown in fig. 2a . At energies just above the K edge tin will attenuate X-rays

more than at energies just below the K edge. The attenuation rises towards lower energies . Thus if a continuous spectrum of X-rays is passed through a tin filter it might be expected that a peak would be produced in the spectrum at energies just lower than the K edge energy of tin, i.e . at energies required for the CT scan of the elements Pd, Ag and Cd . Using known attenuation coefficients for tin and primary spectral shapes for a tungsten target, computer calculations were performed to determine the spectral shape of a pencil beam of X-rays after passing through various thicknesses of tin for various tube voltages. Not only was a more appropriate spectral shape required, the beam intensity had to be reduced to levels which were within the count-rate capabilities of the detection electronics. An attenuation of 102 in the X-ray beam intensity was required . Because the total beam intensity changed more rapidly than spectral shape with increasing filter thickness, it was possible to attenuate the X-ray beam by the required amount and also get a good spectral shape using only tin. Fig. 2b shows the effect of filter thickness on the spectral shape and beam attenuation. The filter chosen was 0.29 mm of tin with an X-ray tube voltage of 40 kV . Also shown is an experimentally measured spectrum for comparison with the calculated spectrum . 3.3 . The X-ray beam collimation and detection count rates A schematic diagram of the apparatus is shown in fig. 3. The X-rays passed through a 1 mm diameter hole in 0.5 mm of lead, then through the tin filter into the scan region . The X-ray source was 18 cm from the lead collimator . The exit aperture was formed from two mutually perpendicular slots in two separate lead plates . The slots were 1 mm wide and the lead 0.5 mm thick.

5

b

N

O

C N

C m â) (1)

o) 3

OY 'L NE

2

O O y lb > vN CY M «Q O

N

0

Sn filters

10

20

30

X-ray energy (keV)

40

U

:r

10

20

30

40

X-ray energy (keV)

Fig. 2. (a) This shows the variation of the mass attenuation coefficient of tin in the vicinity of its K edge at 29 .19 keV. (b) This shows the calculated effect of a tin filter on the X-ray spectrum, as follows (where, for each case, the first number gives the filter thickness ;(c) 0 .3 mm, 4 X 10 -Z ; (d) 0.5 mm, 8.4 X 10 -3. The and the second the calculated transmission) : (a) none, 1; (b) 0 .2 mm, 9X10-2 dash-dotted line is the measured spectrum for 0.29 mm of tin.

J. Fryar et al / Multielement imaging

674

Thus the exit aperture was a square hole of side 1 mm .

200

With a 10 2 attenuation of the beam the count rate was about 10000 per second in the unobstructed X-ray

175

beam with the tube operating at 40 kV and 5 mA . This maximum count rate was determined by the count-rate

150

capability of the multichannel analyser mentioned in section 3.5 .

125

3.4. The scan mechanism The detector used in these experiments was a Si(Li)

detector . Such a detector requires cooling by liquid nitrogen. It is not practical to scan such a device over

N

100

0

75

c ô

the object so the object was moved instead. The scan

50

between 0 ° and 175.5 ° .

25

consisted of 40 X 1 mm steps at 40 rotational positions

X-rays were detected using an EG&G Ortec model

SLP-06180 Si(Li) detector with an energy resolution of

360 eV at 25 keV. The preamplifier of this detector was modified by EG&G so that it could run at count rates

up to and over 10 5 per second . The pulses from the detector were amplified in an Ortec model 575 amplifier to

an

Ortec model 7100

Cd

Pd

Ag y

a

3.5. The detection apparatus

and passed

X-ray absorption spectrum

multichannel

analyser (MCA). At each scan position this MCA accumulated an X-ray spectrum . This MCA had previously

20 21 22 23 24 25 26 27 28 29 30 31 X-ray energy (keV) Fig. 4. This figure shows an example of an X-ray spectrum after transmission through the object . The dips m the upward slope are due to increased absorption near the K edge of the analyte elements Pd, Ag and Cd . The K edges are marked with the small arrows. Shown also are the regions of interest used for the purposes of extrapolation. The regions are about ten channels or 300 eV wide .

been calibrated in energy by using a small radioactive variable energy X-ray source type AMC.2084 made by

Amersham International. With the MCA calibrated in energy it was a simple matter to set up two regions of interest, each about 300 eV and 10 channels wide, on either side of the absorption edge energy for each

element of interest. Eight regions of interest were re-

quired for the three elements . This is because data taken

between absorption edges can be extrapolated to both the low and high energy absorption edge . Fig. 4 shows a spectrum taken at one scan position . The positions of the eight regions of interest are marked . The MCA was the main bottle neck in the data collection . It showed a

percentage dead time of 14% for count rates of 10000 counts/s . The dead time could be improved by using

the low level discriminator of the MCA to remove

counts to the low energy side of the lowest region of

interest . As shown in section 2, corrections for dead-time are not required .

After the spectrum had been accumulated at each

scan position, the total counts in the regions of interest

were read out to a BBC microcomputer. This microcomputer also controlled the scanning of the specimen . At

the end of the complete scan the data was passed to a VAX11/785 computer in which image reconstruction was performed.

4. Image reconstruction Fig. 3. A schematic diagram of the apparatus: (a) 1 mm circular aperture, (b) 1 mm square aperture, (c) object to be scanned, (d) 60 kV X-ray tube, (e) 0.29 mm of tin, (f) Si(Li) detector.

The data stored were counts at each beam position for each of the eight regions of interest . Thus as well as calculating equivalent thicknesses for each of the ana-

67 5

J. Fryar et al. / Multielement imaging

lytes by using eq. (7), it was possible to calculate the quantity In(

N/NO )

= f,u'(r)

dy

(10)

at each region of interest . Here u'(r) is the linear attenuation coefficient at the point r and integration is along the line joining the source to the detector . Therefore it was possible to produce a conventional CT image, showing linear attenuation coefficients, for almost monochromatic X-rays in each region of interest . By using the equivalent thickness scan data of eq . (7) the reconstructed image shows the concentration of the element to whose absorption edge the extrapolation is made. Two methods of reconstruction were used . The first was ART as described in ref. [7]. The algorithm used was, in the notation of ref. [7], C., (new) = C., (old) + Rf,,,e, where e - (t aO - ta

- `fk) 2.

jk

tap and ta; are the equivalent thicknesses as determined from the object and image respectively for the sth linear step on the rth rotation of the specimen . R is a relaxation constant which suppresses oscillations of the concentration with successive iterations . f k is the fractional overlap of the r, s beam with the j, k pixel of the reconstruction grid. This method, although slow to implement on the computer, was the method used to produce quantitative results. The second method used was filtered back-projection using a Ram-Lak filter. Each of the 40 x 40 pixels was divided into quarters and the X-ray beam was divided into two halves (mathematically speaking). If

2 Z

ô Z

c

1

a

Rotation 39 RoI 2

the centre line of a beam half passed through a pixel quarter the filtered projection data was added to that quarter. This produced images in less than 30 s and was used to check the imaging before running the quantitative ART programme which took several minutes for each iteration. This was because of the use of fractional overlap areas to 2% accuracy . 5. Discussion of results The data stored by the microcomputer for passing to the VAXII/785 computer were the counts in each region of interest at each X-ray beam position. As described in section 4 it was possible to reconstruct tomographs showing the variation of the linear attenuation coefficient at each region of interest or the concentration of each analyte. Fig. 5 shows a set of projections for one scan angle. Fig. 5a shows a projection obtained by using eq . (10) on the data for the region of interest 2 in fig. 4. There is no elemental discrimination in this projection and therefore it shows all three analyte elements plus any other elements present. Fig. 5b shows the projection obtained from eq . (1) (but without taking account of the unknown matrix term) using the data in the regions of interest 2 and 3 just to the low and high energy sides of the K edge of palladium respectively . Because no matrix correction has been made the peaks for palladium sit on a negative background . Also shown in fig. 5b are projection data obtained by extrapolation to the K edge of palladium using eq . (7). The negative background has been removed thus showing the usefulness of the extrapolation method . It should be emphasised here that by matrix is meant all elements other than the element to whose K edge the extrapolation is being made ; that includes the other analyte elements. On average the

Y U YN .. E

aa) m Y 7Q

W

10

b

Palladium rotation 39

5

0 -5

1 10 20 30 40 Linear scan position (1 mm steps) Linear scan position (1mm steps) Fig. 5. (a) Shows projection data obtained by use of eq . (10) on the counts in region of interest 2 in fig. 4. These are conventional CT projection data but for a narrow spectral region. From thus type of data tomographs as shown in fig. 6 are produced. (b), (1) A projection showing the variation of the equivalent thickness of palladium. These data were obtained by using eq. (1) on the counts m the regions of interest 2 and 3 of fig. 4. The part of the equation concerned with the matrix has been ignored; hence the negative background . (2) Thus is an equivalent projection to 1 but obtained from eq . (7) which uses extrapolated data. The negative background has been effectively removed.

67 6

J. Fryar et al / Multielement imaging PALLADIUM

_L

Fig. 6. This is an image which shows the variation of the linear attenuation coefficient measured within one region of interest (actually region 2 in fig. 4) . A typical projection used is shown in fig. 5a . Images such as this are useful in showing the overall spatial relationships within the object . The test tubes which held the analytes are clearly seen . There is no elemental discrimination . number of detected photons in each region of interest was about 10 6 per linear scan . From projections as shown in fig. 5 it is possible to reconstruct two types of tomographs. From projections such as shown in fig. 5a are produced tomographs which show the variation of linear attenuation coefficient at an "average" energy within a region of interest . An example of this is shown in fig. 6. The region of interest was that marked 2 in fig. 4. Such an image is useful in seeing the overall spatial relationship between the components of the body . Clearly seen in fig. 6 are the test tubes used to hold the samples. There is some distortion of the two left most test tubes, probably due to the limited number of linear scans. There is some indication of which test tubes hold the higher concentrations of the analytes but there is no indication as to which test tubes hold a particular analyte. From projections such as that shown in fig. 5b are produced tomographs which show the spatial variation of the concentration of a particular analyte. Shown in fig. 7 are the final results of this work ; a set of three images showing the concentration of the three elements Pd, Ag and Cd . The elements are clearly distinguished. Because the mixing, especially of the silver in the form of silver nitrate, was not uniform the concentrations as measured by CT are not necessarily the same as the average value determined from the known quantities of the analytes used . However, the concentrations measured are compatible with the amounts of analyte used . Test tubes containing low and high concentrations are correctly identified . Table 1 contains a comparison of the average concentrations with those obtained in the CT scan .

SILVER

CADMIUM

Fig. 7. A series of images which shows the three elements Pd, Ag and Cd separately . The measured concentrations are com parable to the amounts of analyte added to the epoxy glue (see also table 1) .

For a low density matrix of the size of that used in these experiments the work of Grodzins [8] indicates a detection limit of approximately 10 -3 g/cm3 for a total

J Fryar et al / Multielement imaging Table 1 Comparison of average concentrations of Pd, Ag and Cd with those obtained by CT scanning Element

Concentration [mg/em3 ] Maximum Expected from from image amounts added

Ag Ag

17 12

25 .4 12.7

Cd Cd

15 10

22 11

number of 10 7 photons. This is in agreement with the pixel to pixel variations in fig. 7. It is interesting to note that increasing the total number of photons to 10 9 only reduces the detection limit to about 10 -4 g/cm3 .

677

of this technique in an apparatus which is capable of imaging several elements in a single scan . The effectiveness of the correction method was demonstrated with the clear and separate imaging of the three consecutive elements Pd, Ag and Cd at concentrations within an order of magnitude of the theoretical limits given by the equations of Grodzins [8]. Unfortunately the rate at which data could be accumulated was severely limited by the multichannel analyser used to record the spectra. The use of several single channel analysers and high count-rate capability counters could substantially reduce the data collection time but at the expense of extra complexity in the system . Acknowledgement Thus work was partly funded by the Research and Postgraduate Studies committee of the National Institute for Higher Education, Dublin .

6. Conclusion Conventional X-ray Cr images show the linear attenuation coefficient of the object . This linear attenuation coefficient is determined not only by the chemical composition of the object but also its density. It is not always possible to say whether a change in an image is due to a change in composition or density. The ability to display the concentration of a chosen analyte within the object resolves the problem. By using differential absorption across a K edge this is possible in theory . However, in practice, it is not straightforward . Unless measurements are made very close to the absorption edge, corrections have to be made to the experimental data to allow for matrix absorption effects. In a previous publication [71 we described a method for making this correction . Herein we described the implementation

References

[2] [3] [4]

[61 [7] [81

R.A . Brooks and G. Di Chiro, Phys . Med. Biol . 21 (1976) 390. R.E. Alvarez and A. Macovski, Phys. Med. Biol . 21 (1976) 733. S.J . Riederer and C.A . Mistretta, Med. Phys . 4 (1977) 474. W.J. Glantschmg and A. Holliday, Appl . Opt. 26 (1987) 983 A.C . Thompson, J. Llacer, L. Campbell Finman, E.B . Hughes, J.N . Otis, S. Wilson and H.D . Zeman, Nucl. Instr. and Meth . 222 (1984) 319. Z.H . Cho and K.S . Hong, Nucl. Instr. and Meth. 227 (1984) 385 J. Fryar, K. McCarthy and A. Fenelon, Nucl . Instr. and Meth . A259 (1987) 557. L. Grodzins, Nucl . Instr. and Meth . 206 (1983) 547.