Snapshot molecular imaging using coded energy ... - OSA Publishing

Snapshot molecular imaging using coded energy-sensitive detection Joel A. Greenberg, Kalyani Krishnamurthy, and David Brady∗ Department of Electrical and Computer Engineering, Duke University, Durham, NC, 27708, USA ∗ [email protected]

Abstract: We demonstrate a technique for measuring the range-resolved coherent scatter form factors of different objects from a single snapshot. By illuminating the object with an x-ray pencil beam and placing a coded aperture in front of a linear array of energy-sensitive detector elements, we record the coherently scattered x-rays. This approach yields lateral, range, and momentum transfer resolutions of 1 mm, 5 mm, and 0.2 nm−1 , respectively, which is sufficient for the distinguishing a variety of solids and liquids. These results indicate a path toward real-time volumetric molecular imaging for non-destructive examination in a variety of applications, including medical diagnostics, quality inspection, and security detection. © 2013 Optical Society of America OCIS codes: (340.7430) X-ray coded apertures; (110.7440) X-ray imaging; (110.3200) Inverse scattering; (110.1758) Computational imaging.

References and links 1. J.-P. Schlomka, A. Harding, U. van Stevendaal, M. Grass, and G. L. Harding, “Coherent scatter computed tomography: a novel medical imaging technique,” Proc. SPIE 5030, Medical Imaging 2003: Physics of Medical Imaging pp. 256–265 (2003). URL + http://dx.doi.org/10.1117/12.479949. 2. G. Harding and B. Schreiber, “Coherent X-ray scatter imaging and its applications in biomedical science and industry,” Radiation Phys. Chem. 56(12), 229–245 (1999). URL http://www.sciencedirect.com/science/article/pii/S0969806X99002832. 3. B. Sun, M. Li, F. Zhang, Y. Zhong, N. Kang, W. Lu, and J. Liu, “The performance of a fast testing system for illicit materials detection based on energy-dispersive X-ray diffraction technique,” Microchem. J. 95(2), 293–297 (2010). URL http://www.sciencedirect.com/science/article/pii/S0026265X10000056. 4. D. O’Flynn, C. B. Reid, C. Christodoulou, M. D. Wilson, M. C. Veale, P. Seller, D. Hills, H. Desai, B. Wong, and R. Speller, “Explosive detection using pixellated X-ray diffraction (PixD),” J. Instrumen. 8(03), P03007 (2013). URL http://stacks.iop.org/1748-0221/8/i=03/a=P03007. 5. A. Dicken, K. Rogers, P. Evans, J. Rogers, and J. W. Chan, “The separation of X-ray diffraction patterns for threat detection,” Appl. Radia. Isotopes 68(3), 439–443 (2010). URL http://www.sciencedirect.com/science/article/pii/S0969804309007611. 6. J. Hogan, R. Gonsalves, and A. Krieger, “Fluorescent computer tomography: a model for correction of X-ray absorption,” IEEE Trans. Nuclear Sci. 38(6), 1721–1727 (1991). 7. S. D. M. Jacques, C. K. Egan, M. D. Wilson, M. C. Veale, P. Seller, and R. J. Cernik, “A laboratory system for element specific hyperspectral X-ray imaging,” Analyst 138, 755–759 (2013). URL http://dx.doi.org/10.1039/C2AN36157D. 8. H. Yan, Y. S. Chu, J. Maser, E. Nazaretski, J. Kim, J. J. Lombardo, and W. K. S. Chiu, “Quantitative x-ray phase imaging at the nanoscale by multilayer Laue lenses,” Sci. Rep. 3, 1307 (2013). 9. G. Harding and J. Kosanetzky, “Status and outlook of coherent-x-ray scatter imaging,” J. Opt. Soc. Am. A 4(5), 933–944 (1987). URL http://josaa.osa.org/abstract.cfm?URI=josaa-4-5-933. 10. B. W. King, K. A. Landheer, and P. C. Johns, “X-ray coherent scattering form factors of tissues, water and plastics using energy dispersion,” Phys. Med. Biol. 56(14), 4377 (2011). URL http://stacks.iop.org/0031-9155/56/i=14/a=010.

#196671 - $15.00 USD Received 30 Aug 2013; revised 4 Oct 2013; accepted 9 Oct 2013; published 17 Oct 2013 (C) 2013 OSA 21 October 2013 | Vol. 21, No. 21 | DOI:10.1364/OE.21.025480 | OPTICS EXPRESS 25480

11. J. Delfs and J.-P. Schlomka, “Energy-dispersive coherent scatter computed tomography,” Appl. Phys. Lett. 88(24), 243506–243508 (2006). URL http://link.aip.org/link/?APL/88/243506/1. 12. G. Harding, M. Newton, and J. Kosanetzky, “Energy-dispersive X-ray diffraction tomography,” Phys. Med. Biol. 35(1), 33 (1990). URL http://stacks.iop.org/0031-9155/35/i=1/a=004. 13. A. Dicken, K. Rogers, P. Evans, J. W. Chan, J. Rogers, and S. Godber, “Combined X-ray diffraction and kinetic depth effect imaging,” Opt. Express 19(7), 6406–6413 (2011). URL http://www.opticsexpress.org/abstract.cfm?URI=oe-19-7-6406. 14. K. MacCabe, K. Krishnamurthy, A. Chawla, D. Marks, E. Samei, and D. Brady, “Pencil beam coded aperture x-ray scatter imaging,” Opt. Express 20(15), 16310–16320 (2012). URL http://www.opticsexpress.org/abstract.cfm?URI=oe-20-15-16310. 15. K. P. MacCabe, A. D. Holmgren, M. P. Tornai, and D. J. Brady, “Snapshot 2D tomography via coded aperture x-ray scatter imaging,” Appl. Opt. 52(19), 4582–4589 (2013). URL http://ao.osa.org/abstract.cfm?URI=ao-52-19-4582. 16. R. D. Speller, J. A. Horrocks, and R. J. Lacey, “X-ray scattering signatures for material identification,” Proc. SPIE 2092, Substance Detection Systems pp. 366–377 (1994). URL + http://dx.doi.org/10.1117/12.171256. 17. L. I. Rudin, S. Osher, and E. Fatemi, “Nonlinear total variation based noise removal algorithms,” Physica D: Nonlinear Phenomena 60(1), 259–268 (1992). 18. P. Getreuer, “Rudin-Osher-Fatemi Total Variation Denoising using Split Bregman,” Image Processing On Line 10 (2012). 19. S. Rutishauser, M. Bednarzik, I. Zanette, T. Weitkamp, M. Brner, J. Mohr, and C. David, “Fabrication of two-dimensional hard X-ray diffraction gratings,” Microelectron. Eng. 101(0), 12–16 (2013). URL http://www.sciencedirect.com/science/article/pii/S0167931712004935. 20. B. W. King and P. C. Johns, “Measurement of coherent scattering form factors using an image plate,” Phys. Med. Biol. 54(20), 6437 (2009). URL http://stacks.iop.org/0031-9155/54/i=20/a=C01. 21. J. A. Greenberg, K. Krishnamurthy, M. Lakshmanan, K. MacCabe, S. Wolter, A. Kapadia, and D. Brady, “Coding and sampling for compressive x-ray diffraction tomography,” Proc. SPIE 8858, 885813 (2013). URL http://dx.doi.org/10.1117/12.2027128. 22. P. Seller, S. Bell, R. J. Cernik, C. Christodoulou, C. K. Egan, J. A. Gaskin, S. Jacques, S. Pani, B. D. Ramsey, C. Reid, P. J. Sellin, J. W. Scuffham, R. D. Speller, M. D. Wilson, and M. C. Veale, “Pixellated Cd(Zn)Te high-energy X-ray instrument,” J. Instrumen. 6(12), C12,009 (2011). URL http://stacks.iop.org/1748-0221/6/i=12/a=C12009.

1.

Introduction

The penetrating power of x-rays makes them useful for non-destructive examination of the interior of an object. While traditional transmission-based x-ray imaging systems record only macroscopic object features, more recent applications require one to measure the molecular structure of an object at each location. The ability to perform spatially-resolved material discrimination significantly enhances the performance of medical [1] and industrial radiography [2] as well as contraband [3] and explosives detection systems [4, 5], for example. To this end, a number of approaches have been developed, including fluorescence [6], hyperspectral transmission [7], phase contrast [8], and coherent scatter imaging [9]. Of these techniques, coherent scatter imaging has been the most promising approach because of its capacity to provide material-specific signatures in optically thick materials. In coherent scatter, x-rays scatter elastically from the electron cloud of an atom or molecule. The interference of x-rays scattered from different locations in the material yields a sensitivity to the microscopic structure of the material and produces a unique signature for a variety of materials [9, 10]. The relationship between the object’s microscopic structure and the properties of the scattered x-rays are given by Bragg’s law q≡

E 1 = sin(θ /2), 2d hc

(1)

where q is the momentum transfer, d is the effective lattice spacing of the material, E is the x-ray energy, h is Planck’s constant, c is the speed of light in vacuum, and θ is the angle of deflection of the x-ray. Based on Eq. (1), there are two isomorphic ways that one can measure the scatter intensity for a given value of q: one can either keep E fixed and vary θ (known as #196671 - $15.00 USD Received 30 Aug 2013; revised 4 Oct 2013; accepted 9 Oct 2013; published 17 Oct 2013 (C) 2013 OSA 21 October 2013 | Vol. 21, No. 21 | DOI:10.1364/OE.21.025480 | OPTICS EXPRESS 25481

angle dispersive x-ray diffraction, ADXRD) or keep θ fixed and vary E (energy dispersive xray diffraction, EDXRD). While both approaches have been shown to yield good results, each requires that one substantially filter the x-rays. For example, one must spatially filter the initial x-ray source for both ADXRD and EDXRD. In addition, ADXRD requires a spectrally-filtered source, whereas EDXRD requires one to spatially filter the scattered x-rays. This requirement, coupled with the fact that only a few percent of x-rays are typically scattered coherently from any object voxel, typically necessitates long scan times for large objects. Recently, O’Flynn et al. [4] used a two dimensional array of energy sensitive detector elements to combine ADXRD and EDXRD into a system that rapidly records the coherent scatter signal from a known object location, but this approach requires that the object location be known a priori. In order to build up a volumetric image of an object’s molecular structure, one must unambiguously locate the origin of scatter. Several techniques have been developed to accomplish this task, such as coherent scatter computed tomography (CSCT) [11], selected volume tomography (SVT) [12], and kinetic depth effect x-ray diffraction (KDEXRD) [13], but each requires the acquisition of multiple measurements, which increases the measurement time, implies the need for moving parts, and precludes time-resolved imaging. MacCabe et al. [14, 15] recently demonstrated a snapshot imaging system known as coded aperture x-ray scatter imaging (CAXSI) in which one combines ADXRD with a coded aperture to acquire simultaneously the range and identity of an object. To detect the scattered photons and observe the imposed code pattern, MacCabe et al. employed a two-dimensional detector composed of energy-integrating pixels. In this paper, we demonstrate depth-resolved, snapshot molecular imaging by using a linear array of energy-sensitive detector elements to leverage simultaneously the benefits of ADXRD and EDXRD. By placing a coded aperture in the path of the scattered x-rays, one enhances the throughput of the system without comprising the imaging performance. This system, referred to as coded aperture coherent scatter spectral imaging (CACSSI), allows one to measure quickly both the location and molecular signature of a target object placed in the path of a pencil beam over a range of ∼ 40 cm. In contrast to CAXSI, CACSSI is faster and yields superior resolution. In addition, the use of a linear rather than an area detector array makes the system more compact and demonstrates that the information content in CACSSI exceeds that of CAXSI. The organization of this paper is as follows. Section 2 develops the physical model for CACSSI, analyzes the impact of system parameters on performance, and discusses the reconstruction strategy used to estimate the object. Sections 3 and 4 describe the experimental setup used and results obtained, and Sec. 5 summarizes our conclusions and discusses future directions. 2.

Physical model

As mentioned above, coherent scatter involves the interference of x-rays scattered from centers within an optical coherence length, which is on the order of several nm for a traditional “incoherent” x-ray source. The differential cross section describing the scatter of an x-ray with energy E into a given solid angle is given by dσcoh (E) re2 = 1 + cos(θ )2 f (q). dΩ 2

(2)

Here re is the classical electron radius and f (q) is the square of the coherent scatter form factor, which modifies the Thompson cross section of a free electron by taking into account the effect of nearby scatters. We consider the case of an incident line of x-rays (i.e., a pencil beam) with spectrum Φ(E) illuminating an object with unknown composition and location (see Fig. 1(a)). X-rays at each energy undergo coherent scatter, which causes their paths to deviate by an angle θ that depends #196671 - $15.00 USD Received 30 Aug 2013; revised 4 Oct 2013; accepted 9 Oct 2013; published 17 Oct 2013 (C) 2013 OSA 21 October 2013 | Vol. 21, No. 21 | DOI:10.1364/OE.21.025480 | OPTICS EXPRESS 25482

on q and E through Bragg’s law [Eq. (1)]. By placing energy-sensitive detector pixels along a radial line centered on the primary x-ray beam (we choose here to place them along the x-axis for simplicity), one finds that the number of photons with energy E scattered into a differential solid angle element centered at position x is given as Z zm 2hcqz dσcoh (E) ∆Ω T (E, x, z)t x 1 − δ E− . (3) y(E, x) = dq dz Φ(E) dΩ z x Here the delta function enforces Bragg’s law, we have made the small angle approximation that θ ∼ x/z (which is valid for the energies and detector locations considered here), ∆Ω is the differential solid angle subtended by the detector, and T (E, x, z) = exp[−µ(E, x, z)] and µ are the position- and energy-dependent transmission function and linear attention coefficient of the object, respectively. We note that this model considers coherent scatter only and does not take into account incoherent (Compton) scatter. This approximation is valid for energy-resolved measurements taken at small angles, where the ratio of measured Compton to coherent x-rays is small owing to the peaked nature of coherent scatter in both energy and angle. Multiple scatter, which is also unmodeled but tends to be small and fairly constant, is accounted for approximately through the consideration of a noise term (see Sec. 2.2). energy-sensitive detectors

a) object

Coherent scatter

Coded aperture t(x)

Δx

x Δw zm

x-ray pencil beam

z

b) source

pinhole collimators detector mask

c) 12 mm

object

Fig. 1. (a) Schematic and (b) photograph of the CACSSI setup, including the x-ray source, collimators, coded aperture, object space, and detector. (c) X-ray transmission image of the coded aperture.

The function t(x) in Eq. (3) represents the coded aperture transmission pattern in the plane of the code (i.e., at z = zm ), which is magnified by a factor of z/(z − zm ) upon propagation to the detector plane. This position-dependent magnification disambiguates the angular origin of the scatter and yields range information about the object, thereby allowing us to recover both z and q independently from a single measurement. To gain insight into the role of the code, one can alternatively view the coded aperture as allowing simultaneous measurement of the scatter in two distinct planes. This results in a measurement of the scatter radiance with a standard irradiance detector and, in the case of coherent scatter, yields range information

through Bragg’s law. As discussed in [14], such a mask should be orthogonal in scale so that objects located at different range locations yield optimally-distinguishable code projections. In the following, we choose a periodic binary code t(x) = [1 + sign(sin[ux])]/2 (see Fig. 1(c) for ease of fabrication, where u is the frequency of the code pattern. To understand better the behavior of the system, one can calculate the impulse response by considering the signal generated by a point scatterer described by the position-dependent form factor f (z, q) = fo δ (z − zo )δ (q − qo ) using Eq. (3). We also make the assumption that T (E, x, z) is the same along all paths through the object (i.e., that it is constant and equal to the total attenuation through the object along, To (E)). While this will result in inaccuracies in the reconstruction for highly-structured materials, it is correct for simple objects and greatly simplifies the model. In addition, we note that this approximation will affect only the relative amplitude of f (z, q) and not the locations of the peaks, which is sufficient for the identification of a broad range of material classes [16]. Making these approximations yields the impulse response 2hcqo zo zm zo t x 1− . (4) h(E, x; qo , zo ) ≈ 2 2 3/2 To (E) fo Φ E = x zo (x + zo ) From the form of Eq. (4), one sees that any single pixel yields information regarding the product qo zo . While SVT relies on collimation to focus on a single value of z to avoid this ambiguity, the code used here modulates the scatter intensity across pixels to introduce additional signal structure that depends solely on zo . Thus, each CACSSI measurement consists of a multiplexed signal containing information about the quantities zo and qo zo , from which one can recover uniquely f (z, q). Finally, we model the expected measured signal by convolving the ideal impulse response with the known detector response. We consider each detector pixel to have a spatial extent ∆x and an energy response described by a Gaussian with a half width at half maximum (HWHM) of ∆E. To evaluate the model numerically, we next consider object and measurement vectors f an y with components f j and yi , respectively, and descritize the impulse response to arrive at the form for the forward matrix H with elements " # Z x − xi E − Ei 2 rect Hi j ≈ dE dx h(E, x; q j , z j ) exp −log(2) , (5) ∆E ∆x where rect(x) is equal to unity for |x| < 1/2 and zero otherwise. We typically choose to sample the object space with a voxel size of 5 mm in z and 0.05 nm−1 in q. By considering only measurement data for energies between 20 and 90 keV, the matrix H ∈ RN×Nz Nq is approximately square where N is the number of measurements. In this form, one can solve the linear problem y = H f for f using the algorithm discussed in Sec. 2.2. 2.1.

Resolution analysis

One can estimate the system performance by considering how the choice of system parameters and uncertainties affect the structure of H . We consider first the achievable range resolution. MacCabe et al. [15] showed that, in the absence of noise, the fractional uncertainty in z for the case of a sinusoidal code is given as ∆z/z = z/zm uD[z, q, Φ(E)], where D is the length of the detector that measures coherent scatter. For the case where the angular extend of the scatter exceeds that subtended by a detector with length L, D = L. In contrast, one finds that D = 2hczq(Emax − Emin )/Emax Emin for the case of a point object at f (z, q) and an incident spectrum with total range [Emin , Emax ] when the detector exceeds the scatter. For the latter scenario, which arises in the experimental configuration described in Sec. 3, this yields a fractional spatial #196671 - $15.00 USD Received 30 Aug 2013; revised 4 Oct 2013; accepted 9 Oct 2013; published 17 Oct 2013 (C) 2013 OSA 21 October 2013 | Vol. 21, No. 21 | DOI:10.1364/OE.21.025480 | OPTICS EXPRESS 25484

resolution

∆z Emax Emin = z 2h c q u zm (Emax − Emin )

(6)

which is independent of z, but depends on the properties of the code, object, and incident spectrum. One can understand this result as follows: the determination of object range is commensurate with distinguishing the spatial frequencies of the projected code pattern on the detector plane due to scatter originating at different locations. Thus, one obtains the best results by sampling the maximum number of modulation periods (i.e., using finer mask features, broader spectrum, etc.) for objects placed near the mask (i.e., at locations where the magnification varies most rapidly with range). Similarly, the uncertainty in the estimated momentum transfer ∆q depends on the uncertainty in the scatter angle δ θ and energy of the measured x-rays ∆E. Using Eq. (1), one obtains ∆q = q

"

" =

∆E E ∆z z

2

2

δθ + θ

2 #1/2

∆E δ x + chqz

#1/2 (7)

√ where δ x = ∆x2 + ∆w2 and ∆w is the beam radius. We assume that a detector pixel is located at the optimum location xopt = (2h c qzδ x/∆E)1/2 , which minimizes ∆q. The momentum transfer resolution is therefore dependent on the range resolution, the choice of object, and the energy resolution of the detector. As this approach represents a combination of ADXRD and EDXRD, one can vary how strongly the performance acts like either type of architecture by choosing accordingly the relative uncertainties of the measured x-rays’ scatter angle and energy. For the situation considered experimentally in Sec. 3, we use a detector with a pixel pitch ∆x=0.8 mm and energy resolution ∆E=4 keV. In order to maximize the measured photon flux without impacting the resolution, we choose a beam width ∆w=0.75 mm. Nevertheless, the uncertainty in the recorded photon’s energy dominates the expected momentum transfer resolution. Figure 2 shows the predicted performance based on these chosen parameters as a function of q for three values of z. In general, we find that ∆z decreases with increasing q and decreasing z and is typically on the order of 10 of mm (see Fig. 2(a)). On the other hand, ∆q increases with increasing q and decreasing z and is on the order 0.1 nm−1 (see Fig. 2(b)). We note that we expect these estimates to present an upper bound on the resolution for optically thin objects. This is due to the fact that other mechanisms (e.g., Compton and multiple scatter) can broaden the angular range of modulated scatter present on the detector and thereby improve the spatial resolution. In addition, the performance of nonlinear inversion algorithms can exceed that predicated upon a simple Rayleigh-like criterion and lead to superior results. 2.2.

Reconstruction algorithm

H f + µ b ) as being Poisson distributed with a mean We model the measurements y ∼ Poisson (H value of H f where H is the forward operator discussed in Sec. 2, f is the object that we wish to reconstruct, and µ b is any unmodeled background scatter from the experimental system. Given µ b ), we exploit a maximum a knowledge of H and a realization of the background b ∼ Poisson(µ posteriori (MAP) estimation method with total-variation (TV) regularization to recover f from y . The choice of TV regularization is motivated by the assumption that the underlying image f is piecewise smooth in the z, q space such that the gradient of this image has a small norm.


Dz @mmD

80

z=200 mm z=400 mm z=600 mm

60 40 20 0

1

2

3 4 5 q @1ênmD

b)

0.30

Dq @1ênmD

a)

0.25

6

0.20 0.15 1

2

3 4 q @1ênmD

5

6

Fig. 2. Expected (a) momentum transfer and (b) spatial resolution as a function of q for different values of z (as determined by Eqs. (6) and (7) for u = 1/2 mm−1 , ∆E=4 keV, [Emin ,Emax ]=[20,125] keV, δ x=1 mm and zm =135 mm.

Specifically, we would like to solve the following optimization problem:

bf = arg min − log p y H H ef + µ b + τ ef

(8)

TV

ef

where " e H f +µb = ∑ − log p y H

!

∑ Hk, j fej + µbk j

k

! + yk log ∑ Hk, j fej + µbk

# − log yk !

j

is the Poisson log-likelihood, ef is the total-variation norm of ef and τ is an user-defined TV parameter that balances the log-likelihood term and the penalization term. Since µ b is unknown, b b , which is obtained using the Poisson we replace µ b in the equation above by its estimate µ TV denoising method [17, 18]:

e b TV , e b + λ µ b b = arg min − log p b |µ µ eb µ

where λ is the parameter that dictates the smoothness of the resulting estimate. We solve the optimization problem in Eq. (8) using an iterative, generalized expectation-maximization (GEM) algorithm that consists of the below E- and M-steps at the pth iteration: y bf (p−1) × H T (p−1) bb H bf +µ E-step: z(p) = (9) T H 1 N×1

i h (p)

M-step: bf = arg min − log p z (p) | ef + τ ef (10) ef

TV

for p = 1, 2, · · · . The multiplication and division symbols in Eq. (9) indicate element-wise operations and 1 N×1 is a length N vector of all ones. The M-step in Eq. (10) is the Poisson TV (0) denoising algorithm discussed in [18]. We initialize the algorithm such that bf = H T y and terminate when

(p−1) (p) bf

bf −

2

(p−1)

bf

2

≤ 10−6 . As illustrated in the following section, this algorithm

performs well in practice and preserves fine spectral features while smoothing in regions where the underlying intensity is homogenous.

3.

Experimental methods

To validate the proposed CACSSI scheme, we realize experimentally the configuration shown in Fig. 1(b). We use an x-ray tube (Varian model G1593BI) with a rotating, Tungsten-rhenium anode and a focal spot size of 0.8 mm, which we operate at 125 kV and 250 mA. To produce a pencil beam, we use several stages of collimation: we first collimate the beam it to a 1 mm slit using 2 inch lead bricks placed 100 mm from the focal spot. We then generate the pencil beam by using a pair of 3 mm thick lead sheets with 1.5 mm diameter holes placed at z=350 and 650 mm from the source. This yields a broadband (i.e., spectrally unfiltered) pencil beam with a radius of ∆w=0.75 mm and angular divergence of ∼ 1 mrad. We then place the object on a manual translation stage oriented along z such that it is illuminated by the pencil beam. We place solid samples directly on the translation stage, but use 10 mm diameter Nalgene containers (composed of low density polyethylene) to hold powder and liquid samples. While we observe some scatter from the container around q=0.1 nm−1 , it is small relative to the scatter from the target object. To modulate the scatter, we place a coded aperture downstream from the object (a distance zm from the detector) and out of the path of the primary beam. The coded aperture consists of a comb-like structure where the features (or comb “teeth”) are composed of lead bars with dimensions of 0.8 × 1 × 20 mm (width × thickness × height). We mount these features on an aluminum backing to provide structural support and to enforce a uniform feature frequency of u = 1/2 mm−1 . While smaller feature sizes potentially yield superior imaging performance, these advantages are partially mitigated by the fact that the code must have sufficient thickness to block the x-rays over all energies without limiting its field of view. In addition, fabricating codes with strongly sub-millimeter feature sizes becomes a challenging task [19]. We measure the scattered x-rays using a linear array of 128 energy-sensitive detector pixels (Multix ME-100 Version 1) oriented along the x-direction and located out of the path of the primary beam. The detector pixels are 0.8 × 0.8 mm in extent and composed of CdTe, which is a direct bandgap semiconductor that operates at room temperature. The detector has a nearlyuniform energy resolution of ∆E = 4 keV over the range of 20 - 170 keV, which is binned into 64 energy channels (each with 2.3 keV width) and read out every 20 ms. Because we operate in the low-flux regime (i.e., we record a few thousand counts/s/mm2 ), we find that the detector is dominated by photon noise and is therefore consistent with our Poisson noise model. In a typical experiment, we acquire both a single background (in the absence of an object) and object scatter image, which are used to reconstruct f (q, z) (as discussed in Sec. 2.2). 4.

Results and discussion

We consider first the case of a simple, point-like object. For a 10 mm thick sheet of high-density polyethylene (HDPE) placed a distance z = 252 mm from the detector plane, Fig. 3(a) shows the raw experimental data obtained after 1 s of illumination. As expected, the scatter consists of a curve in (E, x) space given by Bragg’s law that is modulated spatially by the coded aperture. Using the procedure described above in Sec. 2.2, we reconstruct the momentum profile of the object at each range location and show the resulting image in Fig. 3(b). To analyze quantitatively this result, we compute the normalized correlation h bf (·,z) , f (·,z) i

r(z) =

b

f (·,z) , f (·,z)

(11)

between the estimated momentum transfer spectrum at each position and the ‘ground truth’ form factor for the target materials f . We construct f by measuring separately the form factor for each material using standard XRD and assigning it to the known object location, and use the #196671 - $15.00 USD Received 30 Aug 2013; revised 4 Oct 2013; accepted 9 Oct 2013; published 17 Oct 2013 (C) 2013 OSA 21 October 2013 | Vol. 21, No. 21 | DOI:10.1364/OE.21.025480 | OPTICS EXPRESS 25487

a)

b)

c)

d)

Fig. 3. (a) Raw scatter data and (b) the full 2-dimensional reconstruction of the object in q-z space for a 10 mm thick sheet of HDPE placed at z=252 mm. (c) Correlation r as a function of position. (d) Form factors obtained via XRD (blue corresponds to fine q sampling and green corresponds to interpolation at the estimated object’s q sampling) and CACSSI (red). The data was acquired using 500 mAs.

notation a(·,z) to denote all of the q values associated with a particular range location. Plotting this result as a function of z, we find that the reconstructed object is centered at z = 252 mm and displays a full width at half maximum (FWHM) of 11.5 mm (see Fig. 3(c). To evaluate the effective spatial resolution, we note that this measured FWHM corresponds to the real object thickness convolved with a Gaussian blurring kernel of width ∆z ≈5 mm. To evaluate the momentum transfer resolution of the system, we look at f (·,z=252) . Figure 3(d) shows that the estimated form factor is in good agreement with the ground truth when sampled with the same q-spacing as used in the reconstruction. Assuming again a Gaussian blur, we find that ∆q ≈0.2 nm−1 . Thus, the estimated object is very similar to the true object, and the recovered imaging performance is in agreement with that predicted in Eqs. (6) and (7) for objects positioned over a range of approximately 400 mm in z. In addition to materials with discrete, narrow peaks, CACSSI also performs well when imaging amorphous solids and liquids. Figure 4(a) shows the estimated form factors obtained for vials containing different concentrations of methanol (MeOH) and water (H2 O). The shift in peak position to lower q and the disappearance of a bump around q = 0.2 as the concentration of methanol increases is consistent with the XRD measurements (not shown for clarity) and allows one to determine the relative components of liquid mixtures. Similarly, Fig. 4(b) shows that one can also distinguish between water and a 50% water-hydrogen peroxide (H2 O2 ) mixture. Thus, the obtained momentum transfer resolution is sufficient for distinguishing a range of materials. We extend these results by imaging multiple objects placed along the beam path. Figures 5(a)-(c) show the correlation map and estimated form factors for 10 mm thick vials containing graphite and aluminum (Al) powder placed at z = 242 and 267 mm, respectively. Similarly, Figs. 6(a)-(c) shows the estimated form factors for a sequence of three objects along the beam #196671 - $15.00 USD Received 30 Aug 2013; revised 4 Oct 2013; accepted 9 Oct 2013; published 17 Oct 2013 (C) 2013 OSA 21 October 2013 | Vol. 21, No. 21 | DOI:10.1364/OE.21.025480 | OPTICS EXPRESS 25488

a)

b)

Fig. 4. Momentum transfer spectra for a 10 mm diameter vial containing (a) a mixture of MeOH and H2 O at various concentrations and (b) H2 O and H2 O2 (measured one at a time). The data was acquired using 500 mAs for the vial placed at z = 267 mm.

path: water (z=226 mm), Teflon (z=242 mm), and Al powder (z=267 mm). The estimated objects appear at the correct locations in both cases, although the momentum transfer resolution degrades with increasing object complexity. This reduction in image quality stems from a reduction in the signal to noise ratio (SNR) due to attenuation of the primary signal and scatter signals as well as an increase in the noise due to Compton and multiple scatter. In addition, the assumption of uniform attenuation fails for the case of objects with macroscopic structure, although it is straightforward to correct this by measuring directly the attenuation signal [20]. Nevertheless, while it is beyond the scope of this paper, we note that the reconstructed form factors are still sufficiently well-resolved to allow for some degree of material classification [21]. Finally, we look briefly at how the system performs as a function of incident flux by reconstructing a single object using a range of different integration times for a fixed source current. We find that the imaging performance remains nearly uniform down to ∼ 100 mAs before declining smoothly. This trend holds for a variety of solids, plastics, and liquids, although the optimal source flux will depend on the content and organization of the object. The observed degradation of image quality stems from the reduction in modulation contrast due to low SNR at smaller values of mAs. Given that we measure less than ∼ 1% of the total coherent scatter, though, this result suggests that it is possible to realize fast (i.e., few second) volumetric molecular imaging using standard x-ray sources in conjunction with additional detector elements. 5.

Conclusions

In summary, we demonstrate a new approach for performing depth-resolved molecular imaging. By using a coded aperture in conjunction with a linear array of energy-sensitive detector elements, we combine ADXRD and EDXRD to realize a snapshot system that is capable of measuring the coherent scatter form factor of complex objects. While we demonstrate the operating principles here for the measurement of a two-dimensional object (one spatial and one material dimension), this approach can be easily generalized to four-dimensional volumetric molecular imaging by raster scanning the pencil beam throughout the object. Alternatively, one can measure the scatter from a fan beam using a two-dimensional spatial/spectral coded aperture and energy-sensitive detector to image a three-dimensional slice of the object (two spatial and one material dimension). In addition to greatly improving the scan time due to increased parallelism and the collection of a larger fraction of the scattered x-rays, the use of a two-


a)

b)

r

c)

Fig. 5. Object reconstruction results for a 10 mm diameter vial of graphite and Al powder placed at z=242 and z=267 mm, respectively. (a) Correlation along the range dimension and thresholded at 0.75 (red and blue correspond to Al and graphite, respectively). Form factors of (b) graphite and (c) Al obtained via XRD (blue corresponds to fine q sampling and green corresponds to interpolation at the estimated object’s q sampling) and CACSSI (red). The data was acquired using 500 mAs.

a)

b)

c)

Fig. 6. Form factors for a 10 mm diameter vial of (a) H2 O at z=226 mm, (b) Teflon at z=242 mm and (c) Al powder at z=267 mm obtained XRD (blue corresponds to fine q sampling and green corresponds to interpolation at the estimated object’s q sampling) and CACSSI (red). The data was acquired using 500 mAs.


dimensional detector allows one to mitigate the effects of texturing in non-isotropic materials [4]. In order to fully exploit the CACSSI architecture, several potential improvements can be made. As mentioned above, combining scatter and transmission measurements to jointly estimate f (r, r) and µ(E, r) will significantly enhance the performance for complex, opticallythick objects. In addition, using an x-ray source with higher peak voltage will result in a larger incident flux and more x-rays at higher energies, where the scatter is most focused and least affected by attenuation. We can also reduce model error and improve the image quality by modeling Compton scatter and self-attenuation. Finally, the development of high-resolution, multi-dimensional energy-sensitive detectors [22] will lead to superior imaging resolution in shorter times, as the system will be better-conditioned for inversion and provide more information per photon. We therefore believe that the synthesis of coded aperture imaging with energy sensitive detection of coherently-scattered photons will lead to real-time volumetric imaging for application in a range of areas. Acknowledgments We gratefully acknowledge the financial support from the U.S. Department of Homeland Security, Science and Technology Directorate under contract HSHQDC-11-C-00083. We also thank Kenneth MacCabe for helpful discussions and Ehsan Samei for use of his technical resources.
