(p.681) Appendix SOME PROBLEMS OF INTEGRAL GEOMETRY ARISING IN TOMOGRAPHY
(p.681) Appendix SOME PROBLEMS OF INTEGRAL GEOMETRY ARISING IN TOMOGRAPHY
Peter Kuchment
Department of Mathematics, Texas A&M University, College Station, TX 778433368, USA
and
Eric Todd Quinto ^{1}
Department of Mathematics, Tufts University, Medford, MA 02155, USA
This appendix is designed to coordinate the abstract nature of this book with more concrete and practical aspects of the Radon transform.
Beautiful mathematical theorems in integral geometry form the foundation of tomography. The goal of tomography is to recover the interior structure of a nontransparent body using external measurements. In this appendix, we will outline a few tomographic problems that draw on the integral geometry in this book, including support theorems, inversion formulas, range theorems, and microlocal analysis. We will try to show how these theorems are used in tomography. We do not intend to present here any concise survey of tomography, rather a collection of examples chosen according to their connection to the topics considered in the book and to the authors’ interest. The bibliography we provide is also far from being complete, and we refer the reader to books, surveys, and collections of articles [ 14 , 48 – 52 , 67 , 72 , 73 , 80 , 124 , 126 , 138 , 139 , 150 , 161 , 170 ] for details and further references on tomography and related problems of integral geometry.
Xray computed tomography (CT) is the most basic type of tomography; the model is the classical Radon transform that integrates functions over lines. However, the applied mathematics is deep and the properties of the transform provide insight into all parts of integral geometry. In this section we will discuss some of the basic ideas in Xray CT and then talk about various types of limited (p.682) data tomography. Our goal is to show some of the ways integral geometry and microlocal analysis can help one understand these problems.
The goal of Xray CT is to get a picture of the internal structure of an object by Xraying the object from many different directions. Let f : R^{3} → R be the density function of the object; that is, f(x) is the density of the object at the point x ∈ R^{3}. Mathematically, the goal of Xray CT is to recover the density f from these measurements. As Xrays travel on a line L from the Xray source through the object to an Xray detector, they are attenuated by the material on the line they travel along (scatter can be neglected in most Xray CT problems). If the Xrays are monochromatic, then the attenuation at x ∈ L is proportional to the density f(x) (with proportionality constant c > 0). If I(x) is the number of Xray photons in the beam when it arrives at x, then the intensity in a small segment of length Δx is decreased by a multiplicative factor of cf(x)Δx, so
Complete tomographic data is Xray data over all lines. In practice, this means data is collected on a fairly evenly distributed set of lines throughout the object (the concept of complete data can be made precise using sampling theory as in [( 124 ), ( 126 )]). In this case, the commonly used reconstruction algorithm is filtered backprojection, which we will now discuss. We start with some standard notation. Let ω = (cosθ sinθ) be a unit vector on the plane and p be a real number. We denote by L(p, ω) = {x ∈ R^{2}x·ω = p} the line perpendicular to ω and p directed units from the origin. We will also often use ω^{⊥} = (− sin θ, cos θ), the unit vector π/2 units counterclockwise from ω. We should mention that there are two standard ways of parametrizing lines, the one we use in this section (“hyperplane” or nonparametric notation) and the one used in Section ( A.3 ) (“parallel beam” notation). The theorems and ideas in each section are easier to express using the corresponding notation.
The nonparametric Radon transform (( A.2 )) of a function f in the plane can be naturally interpreted as a function of (p, ω) (see Chapter ( 2 )):
(p.683) We define now the backprojection operator, the dual Radon transform, as
Theorem A.1 [( 124 ), ( 151 ), ( 160 )] Let $f\in {C}_{0}^{\infty}\left({R}^{2}\right)$. Then $f=\frac{1}{4\pi}R*\left({I}^{1}Rf\right)\left(x\right).$
The theorem is applied in practice by truncating and smoothing the multiplier σ in I^{−1} and writing this truncated multiplier as a convolution operator in p [( 124 ), ( 126 ), ( 151 )]. The resulting approximate inversion algorithm becomes f ≈ R*(Φ_{*p } Rf) where Φ is the inverse Fourier transform of the truncated multiplier and *_{p} denotes convolution in the p variable.
Proof The proof uses some of the key elementary formulas for the Radon transform. It begins with the projectionslice theorem (which is discussed in Chapter ( 2 ) (( 2.3 ))ff.):
(p.684) To finish the proof, one writes the twodimensional Fourier inversion formula in polar coordinates and uses ( A.3 ) to get
In Section ( A.2.2 ) we will use the properties of R as an elliptic Fourier integral operator (FIO) to develop subtle results about singularity detection in tomography. Now, we will use the projectionslice theorem ( A.3 ) to show R is, in fact, an elliptic FIO. Let J :C ^{∞}(R ^{2}) → C ^{∞}(R×S ^{1}) be defined by Jf(p,ω) = f(pω). We take the onedimensional inverse Fourier transform in ( A.3 ):
A similar exercise using the fact that R* Rf(x) = (2/x)* f shows that
Wave front sets provide a subtle and elegant classification of singularities of functions and distributions, and microlocal analysis allows us to understand how FIOs transform wave front sets (compare Section ( 5.3 )).
In particular, we will use this methodology to learn about how the FIO R detects singularities and to predict how well reconstruction algorithms will perform when they have limited tomographic data. The wave front set of f ∈ 𝒟′(R^{ n }) is defined in Chapter 5 . It is useful because it classifies not only singularities of f using points x, _{0} at which f is not smooth, but also directions above each point in which f is not smooth. We also need the concept of Sobolev wave front set. The distribution f is in H ^{9}(R^{2}) if and only if its Fourier transform, f̃ = ℱ_{ x→ξ} f, (p.685) is in L^{2}(R^{2}, (1 +ξ^{2})^{ s }). This relates global smoothness of f to integrability of its Fourier transform. A local version of this at a point x _{0} ∈ R^{ n } is obtained by multiplying f by a smooth cutoff function $\mathrm{\psi}\in {C}_{0}^{\infty}\left({\text{R}}^{n}\right)$ (with ψ(x _{0}) ≠ 0) and seeing if $\tilde{(\psi f)}$ is in this weighted L ^{2} space. However, this localized Fourier transform $\tilde{(\psi f)}$ gives even more specific information—microlocal information— namely, the directions near which $\tilde{(\psi f)}$ is in L^{2}(R^{2},(1 +ξ^{2})^{ s }). The precise definition is (see [ 135 ], p. 259)
Definition A.2 A distribution f is in the Sobolev space H ^{ s } locally near x _{0} ∈ R^{ n } if and only if there is a cutoff function ψ ∈ D(R^{ n }) with ψ(x_{0}) ≠ 0 such that the Fourier transform $\tilde{(\psi f)}$(ξ) ∈ L^{2}(R^{ n },(1 +ξ^{2})^{ s }). Let ξ_{0} ∈ R^{ n } \ 0. The distribution f is in H ^{ s } microlocally near (x_{0},ξ_{0}) if and only if there is a cut off function ψ ∈ D(R^{ n }) with ψ(x_{0}) ≠ 0 and function u(ξ) homogeneous of degree 0 zero and smooth on R^{ n } \ 0 and with u(ξ_{0}) ≠ 0 such that the product $u\left(\xi \right)\tilde{\left(\mathrm{\psi f}\right)}\left(\xi \right)\in {L}^{2}({R}^{n},(1+\xi {}^{2}{)}^{s})$. The H ^{ s } wave front set of f, WF^{ S }(f), is the complement of the set of (x_{0}, ξ_{0}) near which f is microlocally in H ^{ s }.
It follows from this definition that, if (x_{0}, ξ_{0}) ∉ WF(f), then for all s, f is H ^{ s } near (x_{0}, ξ_{0}). The wave front set and microlocal Sobolev smoothness are usually defined on T*(R^{ n }) \ 0, the cotangent space of R ^{ n } with its zero section removed because such a definition can be extended invariantly to manifolds using local coordinates. To this end, let x_{0} ∈ R^{2}. If y = (y_{1}, y _{2}) ∈ R^{2}, then we let ydx = y _{1}dx_{1} +y_{2}/_{2}dx_{2} be the cotangent vector corresponding to y in ${T}_{{x}_{0}}^{*}{R}^{2}$. Now let (p,ω) ∈ R × S^{ l }. We let dp and dω be the standard basis of ${T}_{(p,\omega )}^{*}(R\times {S}^{1})$. Here, θ maps to ω = (cosθ, sinθ) provides coordinates on S ^{1} and the basis covector dω is the dual covector to ∂/∂θ. The wave front set is extended to distributions on R × S^{ l } using these local coordinates, and it is a subset of T*(R × S^{1}).
The fundamental theorem that gives the relation between Sobolev wave front of a function and its Radon transform is the following.
Theorem A.3 [( 146 )] Let $f\in {\epsilon}^{\prime}\left({R}^{2}\right),{x}_{0}\in L({p}_{0},{\omega}_{0}),{\eta}_{0}=\mathrm{dp}({x}_{0}\xb7{\omega}_{0}^{\perp})\mathit{d\omega},$ and $a\in {R}^{n}\backslash 0.$ The correspondence between WF(f) and WF(Rf) is
This theorem is a corollary of fundamental results of Guilleman [( 59 )], [( 61 )] on the microlocal analysis of R. The general results will first be described and then the relation to Theorem A.3 will be given and finally the proof will be outlined.
Guillemin first developed the microlocal analysis of the Radon transform. For general Radon transforms defined by double fibrations on manifolds with smooth (p.686) nowherezero measures, he proved that R is an elliptic FIO. He also proved that R*R is an elliptic pseudodifferential operator under the Bolker assumption on the Lagrangian manifold associated to the FIO R [( 60 )]. At that time, he developed a new description of FIO using push forward and pull backs generalizing those used to define Radon transforms using double fibrations [( 59 )]. This work was included in [( 61 )] (see pp 364ff.).
Quinto [( 140 )] described the symbol of R*R for all generalized Radon transforms satisfying the Bolker assumption in terms of the measures involved into the transform. He also proved more results for the hyperplane transform, including a simple expression for the symbol and a proof of invertibility for translation invariant transforms. Beylkin [( 19 )] considered Radon transforms satisfying the Bolker assumption integrating over surfaces in R ^{ n } and analyzed the asymptotics of the symbol of R†KR where R† is a backprojection related to Rη and K is a type of filter. We will discuss the generalized Radon transform further in Section ( A.2.3.3 ) and Section ( A.3. )
The microlocal correspondence between WF(f) and WF(Rf) in Theorem A.3 follows from Guillemin's seminal work [( 60 ), ( 61 )]. Because R is an FIO corresponding to a specific Lagrangian manifold, ( A.6 ) holds. Sobolev continuity is a basic property for FIOs; any FIO, such as R, of order ½ will map functions in H ^{ s } of fixed compact support continuously to functions in ${H}_{\mathrm{loc}}^{s+1/2}$. If the operator is elliptic, then the function must be ½ order less smooth than its image. Theorem A.3 is a refinement of this observation. Not only does R smooth of order ½ but it maps functions that are in H ^{ s } microlocally near a given covector to functions that are in H ^{s+½} near the covector given by the correspondence ( A.6 ).
The proof outline is as follows. For this classical Radon transform, the key to the correspondence ( A.6 ) is ( A.4 ), that R is an elliptic FIO with phase function φ(x, p, ω, ω) = (p − x · ω)ω and the explicit Lagrangian manifold associated to this operator [( 60 ), ( 61 )] (see [( 140 )] for details). The fact that (p _{0}, ω_{0}; aη_{0}) in ( A.6 ) determines (x _{0}, aω_{0}dx) is a simple exercise. The assertion about H ^{ s } is given in [( 146 )], and, as discussed above, is proven using general Sobolev continuity results for FIO [ 76 , 168 ].
Theorem A.3 allows one to understand what R does to singularities in a precise and rigorous way and it provides an application of microlocal analysis to tomography. It gives an exact correspondence between singularities of f and those of Rf. Moreover, it states that the singularities of Rf that are detected from the data are of Sobolev order ½ smoother than the corresponding singularities of f. For typical singularities of f (jump singularities in H ^{½}∈) one can realistically expect the corresponding singularities of Rf not to be masked by noise.
The theorem has the corollary
Corollary A.4 Data Rf for (p, ω) arbitrarily near (p _{0}, ω_{0}) detects singularities of f perpendicular to the line L(p _{0}, ω_{0}) but not in other directions. If Rf is in H ^{s+½} locally near (p _{0}), then f is microlocally in H ^{ s } at all points on (p.687) L(p _{0}, ω_{0}) at all directions conormal to L(p _{0}, ω_{0}). Conversely if a direction above (p _{0}, ω_{0}) is in WF^{s+½}(Rf) then f has H ^{ s } wave front set conormal to L(p _{0}, ω_{0}) at the point given by ( A.6 ).
Corollary A.4 follows from the correspondence ( A.6 ) for WF and WF^{ S } and the fact that data WF^{s+½}Rf above (p _{0}, ω_{0}) provides information only about singularities of f conormal to L(p _{0}, ω_{0}). Palamodov stated a closely related idea in [( 133 )]. The “tangent casting” effects of [( 156 )] are related to Corollary A.4.
A simple illustration will give an intuitive feeling for the corollary. Let f : R ^{2} → R be equal to one on the unit disk and zero outside. Then, $Rf(p,\omega )=2\sqrt{1\mathrm{}{p}^{2}}$ for p ≤ 1 and Rf(p, ω) = 0 for p > 1. The only lines where Rf is not smooth are those with p = 1 and these lines are tangent lines to the unit circle. Since WF(f) is the set of vectors normal to the unit circle, singularities of Rf precisely capture singularities of f.
Limited tomographic data is tomographic data given on some proper open subset of R × S ^{1}. Theorem A.3 and Corollary A.4 provide a paradigm to decide which singularities of f are stably visible from limited tomographic data, and we will examine what this predicts for the three common types of limited data: limited angle CT (Section ( A.2.3.1 )), exterior CT (Section ( A.2.3.2 )), and the interior problem or local CT (Section ( A.2.3.3 )).
The basic idea is as follows. If Rf is in H ^{s+½} near (p _{0}, ω_{0}) ∈ R × S ^{1}, then, by Corollary A.4, f must be in H ^{ s } at all points of L(p _{0}, ω_{0}) in directions conormal to L(p _{0}, ω_{0}). By Theorem A.3, if Rf is not in H ^{s+½} microlocally in some direction at (p _{0}, ω_{0}) then f cannot be in H ^{ s } conormal to L(p_{0}, ω_{0}) at a specific point on this line that is given by ( A.6 ). So, singularities of Rf above (p_{0}, ω_{0}) determine specific singularities of f conormal to L(p _{0}, ω_{0}). However, data Rf arbitrarily close to (p_{0}, ω_{0}) does not tell anything in a stable way about singularities of f away from L(p _{0}, ω_{0}) or at points of L(p _{0}, ω_{0}) in directions not conormal to this line. By ( A.6 ), these singularities are visible from other data.
A typical tomographic density f is often a piecewise smooth function that is smooth on open sets with wellbehaved boundaries. So, singularities of f occur at the boundaries, and the singularities are in H ^{½ − ∈} for ∈ > 0. By Theorem A.3, singularities of the Radon data Rf will be in H ^{l − ∈}. A limitation of this analysis is that any discrete data Rf can be considered smooth. However, standard singularities of Rf are going to have large norm in H ^{1} and so should be visible. Moreover, the paradigm of the preceding paragraph is observed in all typical CT reconstructions from limited data.
Ramm and Zaslavsky [( 152 )] have analyzed how the Radon transform itself behaves on functions that are smooth except at smooth boundary surfaces. They give precise asymptotics of Rf at lines tangent to boundaries depending on the curvature of the boundaries, and they have proposed a singularity detection (p.688) method using this information on the raw data [( 152 )]. This method has been tested on simulated data [( 84 )]. A more general method has been proposed [( 146 )].
We now apply Theorem A.3 and Corollary A.4 as well as some classical theorems to three common types of limited data tomography problems in the plane. We will use the microlocal results to determine which singularities of f should be visible from limited tomographic data. One can also understand stability of limited data tomography using singular value decompositions [( 35 ), ( 107 ), ( 110 ), ( 113 ), ( 114 ), ( 116 )], and they reflect the principle predicted by the microlocal analysis: singular functions associated with large singular values (which are easy to reconstruct) oscillate in directions in which the wave front is easily detectable, and vice versa.
A.2.3.1 Limited angle tomography
Let U be a proper open subset of S ^{1} such that U = −U. Limited angle data is tomographic data given on R × U only, and the limited angle Radon transform is the Radon transform $R:{L}_{0}^{1}\left({\mathbf{\text{R}}}^{2}\right)\to {L}_{0}^{1}(\mathbf{\text{R}}\times U)\mathrm{where}{L}_{0}^{1}$ is the set of L ^{ l } functions of compact support. The goal is to recover f using data Rf(p, ω) for (p, ω) ∈ R × U. This problem comes up in electron microscopy [( 33 )] and in Xray airport inspections of luggage traveling on a belt. In both cases, one cannot acquire Xray CT data all around the object but only in a limited angular range. Alternatively, in electron microscopy, one has a picture that shows many molecules with the same crystalline structure but oriented in random directions. This is like regular CT with complete data, but with no information about the directions [( 126 )].
Using the projectionslice theorem ( A.3 ), one can prove this limited data transform is injective for functions of compact support: knowing Rf on R × U gives the Fourier transform f̃ on the cone in R ^{2}, V = {tωt ∈ R, ω ∈ U}. Since f has compact support, f̃ is real analytic and this Fourier data determines f. Tuy [ 169 ] developed a reconstruction algorithm for this problem by extending this Fourier data f̃ on the cone V to R ^{2} by making f̃ equal to zero off of V, then taking the inverse Fourier transform. This method destroys compact support since the Fourier transforms of compactly supported functions are real analytic. His inversion is really a pseudodifferential operator that is elliptic in directions from V so standard results about wave front sets [( 76 ), ( 168 )] show that all wave fronts of f in the directions in U will be preserved.
Other appealing and successful inversion methods involve extending the Radon data from R × U to R × S ^{1} using the range conditions on the Radon transform [( 103 ), ( 107 )]. The basic idea is to project the data onto the range of the Radon transform with complete data (on R × S^{1}) and then invert this consistently completed data using filtered backprojection. Here are the details. Let B be a ball of radius M > 0. To implement this, one uses the singular value decomposition for R : L ^{2} (B) → L ^{2} ([− M, M] × S ^{1}) in appropriate measures [( 30 ), ( 124 )]. One projects the incomplete data onto this basis of singular functions on L ^{2} ([− M, M] × S ^{1}). (p.689) The projection is unstable because it requires an extrapolation: the singular basis is not orthogonal on [− M, M] × U but on the larger set [− M, M] θ S ^{1}. The condition number of the projection matrix increases exponentially as the set Ubecomes smaller [( 35 ), ( 107 )]. Alternatively, if one were to construct an orthonormal basis on L ^{2} [− M, M] × U then the singular functions would have large norm outside of U.
This instability is intrinsic to the problem, and it is reflected in the fact that some singularities of f are not stably visible from the limited angle data. By Corollary A.4, the only singularities of f that can be detected in a stable way are those with directions in U. To see this, choose x ∈ R ^{2} and ω ∈ U. Any wave front of f at (x; ωdx) is detected by limited angle data because the line L(x·ω, ω) is in this data set. For the same reason, the wave front of f at (x; ωdx) for ω U will not be stably detected by this limited angle data.
Moreover, an argument of Finch related to this singularity analysis shows that inversion of limited angle or exterior data is illposed in any range of Sobolev norms [( 40 )]. The argument for limited angle CT is as follows. We construct a nonsmooth function f that is smooth at all points in R ^{2} in all directions in U. Then, by ( A.6 ), Rf is smooth on R × U. So, Rf is in every local Sobolev space above R × U but f is not in all Sobolev spaces on R ^{2} since f is not smooth. Therefore, inversion of the limited angle Radon transform is not continuous in any range of Sobolev norms.
Here is one way to construct such an f. Let ω ∉ U, let L be the line through the origin normal to ω and let H be a halfplane with boundary L. We take a smooth nonzero radial function of compact support, ψ and let f be the result of setting ψ to zero on H. Then, f is smooth at all points in R ^{2} in all directions in U but f is not in C ^{∞}. In fact, WF(f) is a subset of points in L with direction ωdω, i.e. conormals to L.
This instability is also illustrated by the singular functions for the limited angle problem [( 107 )]. Those corresponding to large singular values (easy to reconstruct) oscillate generally in directions in U and those corresponding to small singular values (hard to reconstruct) oscillate generally in directions outside of U.
A.2.3.2 The exterior problem
Let M > 1 and assume supp f ⊂ {x ∈ R ^{2}x ≤ M}. In the exterior problem, one has data Rf(p, ω) for all ω but only for p > 1. By the support theorem for the Radon transform ([( 67 )], see Chapter 2 ), one can reconstruct f(x) for x > 1. This problem comes up in studies of the solar corona [( 6 )] in which data is total intensities of the corona of the Sun along lines exterior to the core of the Sun that go from the solar corona to the observer on Earth. Exterior data also occurs in industrial tomography of very large objects, for which Xray data through their centres is too highly attenuated to be usable [( 148 )].
Let x > 1 and ω ∈ S ^{ l }. Then the only singularities of f at x that are reconstructed in a stable manner are those for ω with L(x · ω, ω) in the data set, i.e. forx · ω > 1. Other singularities of f are not stably detected. Because (p.690) some wave front directions are not stably detectable by exterior data, inversion for the exterior problem is highly illposed; a translate of Finch's example in Section ( A.2.3.1 ) can be used to show inversion of the exterior transform is discontinuous in any range of Sobolev norms. Lissianoi [( 105 )] has extended Finch's Sobolev discontinuity result to show that in the exterior problem even recovery of the function in a smaller ring than where the data is given does not help to improve stability. However, logarithmic stability has been proven by Isakov [( 78 )].
Lewitt and Bates [( 103 )], Louis [( 108 )], and Natterer [( 120 )] have developed good reconstruction algorithms that use exterior data. Lewitt and Bates’ algorithm completes the exterior data by projecting it on the range of the complete Radon transform as discussed in Section ( A.2.3.1. ) The projection step is unstable because the singular functions are not orthogonal on the annulus but on a disk. Natterer's algorithm is an effective regularization method. Quinto has developed an exterior reconstruction algorithm which employs a singular value decomposition [ 134 ] for the Radon transform on domain L ^{2}({x ∈ R ^{2}x ≥ 1}) and a priori information about the shape of the object to be reconstructed. Reconstructions for “medical” phantoms are in [( 142 )] and those for industrial phantoms are in [( 143 ), ( 144 )] and real industrial data in [( 148 )]. Exactly those singularities that are supposed to be stably reconstructed are clearly imaged in the reconstruction. In Quinto's reconstruction algorithm, singularities that are not “visible” are blurred.
A.2.3.3 The interior problem and local tomography
Let M > 1 and assume supp f ⊂ {x ∈ R ^{2}x ≤ M}. Interior tomographic data is data Rf(p, ω) for all ω but only for p < 1. Data is missing over lines outside the unit sphere, even though supp f can meet the annulus {x ∈ R ^{2}  l ≤ x ≤ M}. The goal of interior CT is to reconstruct information about f(x) for x ≤ 1. This problem comes up whenever scientists want information only about some portion of the object to be reconstructed, not the whole object, or in problems, such as highresolution tomography of very small parts of objects or electron microscopy, for which it is difficult or impossible to get complete highresolution CT data.
Simple examples (derived using the range theorem for the Radon transform) show the interior transform, the Radon transform with this limited data, is not injective. However, according to Corollary A.4, one can detect all singularities of f in x < 1. To see this, choose a point x inside the unit disk and choose a direction ω ∈ S ^{1}. Then the line through x and normal to ω is in the data set for interior tomography. Therefore, by Theorem A.3, any singularity of f at (x; ωdx) is stably detected by interior data.
Maaβ [( 114 )] has developed a singular value decomposition for this problem. See also [( 110 )]. The authors of [( 103 )] have developed a reconstruction method that projects the interior data onto the range of the Radon transform with complete data and then inverts this completed data. Because of the nonuniqueness of the problem, this projection step is not unique. However, Maaβ has shown that singular functions associated to small singular values are fairly constant inside (p.691) the region of interest, the unit ball [( 114 )]. These singular functions corresponding to small singular values are difficult to reconstruct, but they do not add much detail inside the region since they are relatively constant there. This reflects the fact that all Sobolev singularities inside the region of interest are stably reconstructed.
Even though the interior transform is not invertible, we have shown that all singularities in the region can be stably reconstructed. This explains why singularity detection methods work so well for this problem. Lambda tomography [( 38 ), ( 39 ), ( 160 ), ( 171 )] is a welldeveloped algorithm that finds singularities of a function using interior data over a very small region. One starts with the filtered backprojection inversion formula, f = 1/4πR* I ^{−1} Rf, and replaces I ^{−1} by I ^{−2} = I ^{−1} ο I ^{−l} = −d ^{2} / dp ^{2}. The beauty of this idea is that the formula
For the Xray transform with sources on a curve in R ^{3}, Louis and Maafi [( 109 )] have developed a very promising generalization of lambda CT. Greenleaf and Uhlmann [( 56 ), ( 57 )] completely analyzed the microlocal properties of R* R for admissible transforms on geodesies, including this example. The microlocal analysis of this threedimensional lambda CT operator has been investigated by Finch [( 42 )], Katsevich [( 83 )], and Lan [( 102 )]. The operator R*λ adds singularities to f because R is not well enough behaved; for the classical line transform in the plane, R*λ R is an elliptic pseudodifferential operator and therefore preserves singularities. The theorem corresponding to Theorem A.3 for the Xray transform on lines through a curve is given in [( 146 )], and it explains which singularities are stably reconstructed by this transform. The prediction is observed in the reconstructions in [( 109 )]. The general analysis in [( 62 )] and [( 56 )] provides the microlocal results needed to prove this theorem. The microlocal properties of this line transform were also studied in [( 22 )] as a basis for support and uniqueness theorems for the Xray transform.
Limited angle and exterior versions of lambda CT have been developed, and they are promising on simulations [( 92 )] and tests on industrial and electron microscope data (E.T. Quinto unpublished).
(p.692) Another way to deal with interior data is to truncate the data to zero by multiplying by a cutoff function which is equal to one in the region of interest. If one now applies a standard Radon inversion formula to the truncated data, one can check that the result is a pseudodifferential operator with the symbol equal to one in the local region of interest. This means that although the exact values of the function are not recoverable from the local data, its singularities are, by simple usage of truncated data in inversion formulas. One sees, for instance, that the values of jumps of the function are recovered correctly [( 86 ), ( 90 )]. This idea was further developed by Katsevich and Ramm [( 85 )] under the name of pseudolocal tomography. The authors of [( 38 )] have provided a different method to calculate the jumps of f at boundaries using lambda CT.
Madych [( 115 )] used wavelet analysis to show the strong relationship between lambda CT and pseudolocal CT and regular filtered backprojection. Authors of the papers [( 36 ), ( 153 )] have also used wavelet techniques for local tomography. They use Radon data to calculate wavelet coefficients of the density to be reconstructed. Although they need some data slightly outside the region of interest, their methods are fairly local.
Candès and Donoho [( 27 )] have defined ridgelets, wavelets that are not radially symmetric and are more sensitive to singularities in specified directions. They have developed a tomographic reconstruction algorithm using ridgelets that provides highquality reconstructions in practice. Their ideas reflect the singularity detection predictions of Theorem A.3.
Technicians currently use the sinogram, the graph of Rf(p, ω) in rectangular coordinates, to find boundaries, but this method is subjective, and some industrial data is too homogeneous for this to work [( 149 )].
Local tomography has been developed for generalized Radon transforms [( 92 )], such as the attenuated Radon transform discussed in Section ( A.3 ) and in particular in Section ( A.3.1. ) If ν = ν(x,ω) is a nowherezero smooth weight, we define the generalized Radon transform of $f\in {C}_{0}^{\infty}\left({R}^{2}\right)$ to be
An important generalized Radon transform is the attenuated transform T _{μ} defined in Section ( A.3 ), ( A.11 ). The weight for this transform is
The generalized dual transform of R _{ν} is
We use this to develop a generalized lambda tomography. The calculus of FIOs shows that the generalized lambda operator (compare with ( A.7 ))
In general, the attenuated Radon transform will have a measure that is not smooth on boundaries of regions at which the attenuation changes. This gives rise to pseudodifferential operators ${R}_{\nu}^{*}{d}^{2}/{\mathrm{dp}}^{2}\mathrm{R\nu}$ that have nonsmooth symbols. Katsevich studied these operators, and he proved that one can detect jumps of f that are not masked by singularities of the measure using this information, if one knows the measure [( 82 )]. This work is important even if one does not know the measure ν exactly. Because the operator is an elliptic pseudodifferential operator (albeit with a nonsmooth symbol), the singularities of f that are not masked by singularities of ν will still be visible from the data, even if the measure is unknown. One practical issue would be whether the large amount of noise in SPECT data might create problems for any singularity detection algorithm that “overdifferentiates” the data.
Some problems of medical tomography, radiation therapy, and industrial nondestructive testing naturally lead to consideration of a special type of weighted Xray transforms. We will briefly describe one such model, referring the reader to publications [( 24 ), ( 31 ), ( 32 ), ( 95 ), ( 124 ), ( 138 ), ( 139 ), ( 158 )] for further details, examples, and references.
The socalled single photon emission computed tomography (SPECT) deals with the situation when a nontransparent body contains a distributed radiation source. The goal is to reconstruct the interior intensity distribution f(x) of the (p.694) radiation sources by conducting exterior measurements of the intensity of outgoing radiation. Let us denote by μ(x) the linear attenuation coefficient (or absorption) of the body at the point x. This means that due to the absorption a beam passing through the point x suffers at a small distance Δx the relative intensity loss equal to μ(x)Δx. Then assuming that effects of scatter are small and hence can be neglected (which is in fact not always true), the stationary onevelocity transport equation can be written as
The operator T _{μ} is said to be the attenuated Xray transform of the function f(x) and μ(x) is called the attenuation. The formula above is a nonparametric version of this transform. One can also introduce a parametric one. As mentioned in the last section, we will use a slightly different parametrization of lines in this section than we did in Section ( A.2 ). The theorems and ideas in each section are easier to express using the respective notation. Namely, let ω be the unit vector parallel to the direction of the oriented line L, ω^{⊥} be its orthogonal hyperplane, and y = L∩ω^{⊥}. Then the oriented lines L are in onetoone correspondence with pairs (ω, y), where y ∈ ω^{⊥}. (Now, this notation differs from that of the previous chapters; in Chapter 6 we studied attenuations which do not depend on y.) In the present notation we write
One notices that in contrast with the standard Xray transform f → ∫_{L} f(x)dx, ( A.2 ), the resulting function depends on the orientation of the line L. (p.695) In practice one often averages over the two orientations, thus arriving at a function of nonoriented lines.
In many algorithms one assumes that the attenuation is constant inside the body and zero outside. Let us also assume that the body has a known convex shape. In this case, according to [( 117 )], the attenuated transform can be reduced to a simpler one. Indeed, the integral
Among the natural questions to ask about the attenuated and exponential Xray transforms are:
Injectivity. Can a function of a natural class be reconstructed from its attenuated or exponential transforms? In other words, are these operators injective?
Inversion formulas in cases when injectivity is established.
Stability of inversion.
Range. Are these operators surjective in natural functional classes? Judging by the precedent of the standard Radon transform, the reader certainly should feel that the answer is probably negative, and some nontrivial range conditions should be satisfied. Then the question arises of describing the ranges of these operators. (p.696) Simultaneous reconstruction of the sources density f and attenuation μ. In most cases not only the value of the intensity f(x) distribution, but the attenuation coefficient μ(x) as well is unknown. Is it possible to extract any information about both functions from the values of T _{μ} f or R _{xμ} f?
These problems will be briefly addressed below. We would also like to mention papers [( 37 ), ( 95 ), ( 131 ), ( 132 ), ( 157 )], where some interesting relations between the properties of the exponential transform and complex analysis are discovered that are beyond the scope of this text.
Uniqueness of reconstruction is obviously the first question one should ask. Nonuniqueness would mean that one is unable to recover the interior structure of the object of interest from the exterior information. When we discuss uniqueness we will assume that the function f to be reconstructed has compact support.
The problem of uniqueness for the attenuated Xray transform in dimension 2 has proven to be hard. The first results on uniqueness were the local ones. Since for smooth μ the operator T _{μ} (after averaging over opposite directions) is an elliptic FIO [( 19 ), ( 60 ), ( 61 ), ( 140 )], one concludes that for each point x there is a neighborhood U of x such that uniqueness holds for functions f with support in U. The idea is that locally the attenuated transform looks almost like the standard Xray transform. It was shown in [( 118 )] that C^{2} smoothness of the attenuation μ is sufficient. The paper [( 66 )] contains the proof of the following statement. Let us require that the supports of all functions f under consideration belong to a fixed bounded domain Ω Assume, for instance, that Ω is the unit ball centreed at the origin. Then the transformation T _{μ} is a semiFredholm operator between the Sobolev spaces ${H}_{0}^{s}\left(\mathrm{\Omega}\right)$ and H ^{s+½}(S ^{1} × [−1, 1]) for any s (see also [( 124 )]). This means that the range of this operator is closed, and its kernel Ker T _{μ} is finite dimensional. In the class of such operators property Ker A = 0 is stable under perturbations of small operator norm. This observation immediately leads to a local uniqueness result: for small Ω, operator T_{μ} is norm close to the Xray transform, and hence has zero kernel. The next step was made by Finch [( 41 )], who managed to obtain a “semilocal” uniqueness results. Namely, he proved uniqueness under the condition μ_{∞} diamΩ < 5.37. This result is sufficient for many practical situations. For instance, in medical applications it restricts the diameter of an object to 35.8 cm. The proof is nontrivial and involves energy estimates. A breakthrough came recently in brilliant works [( 8 ), ( 128 ), ( 129 )], where a positive solution of the uniqueness problem under some mild smoothness condition on the attenuation was obtained [( 8 )] and an explicit inversion formula found [( 128 ), ( 129 )] (see also [( 26 ), ( 43 ), ( 58 ), ( 96 ), ( 98 ), ( 125 )] for different derivations and implementations). We will discuss the inversion formulas in a little more detail in the next section.
It was noticed in [( 41 )] that a similar uniqueness problem in dimensions 3 and higher is trivial. The reason is that the threedimensional attenuated Xray (p.697) transform is the collection of the corresponding twodimensional transforms in all affine planes. If there is a compactly supported function f(x) annihilated by the transform, then we can choose a plane which just barely touches the support of f. In this plane one can use the local uniqueness theorem to conclude that the function is in fact zero there. Proceeding further in this way, one can “eat the support away” and finally conclude that it is empty.
The problem of uniqueness was also considered for transforms with more general positive weights ν the generalized Radon transform (A.8). Here there are local uniqueness results for fairly arbitrary weights and global uniqueness theorems for some special types of weights (translation invariant, rotation invariant) [( 100 ), ( 140 ), ( 141 )]. In the case of an analytic weight uniqueness in the class of compactly supported functions follows from analytic ellipticity of the corresponding FIO [( 21 )]. Analytic ellipticity has been used to prove support theorems and uniqueness in many settings such as for Radon transforms on geodesic spheres in realanalytic manifolds [( 145 )] and for circular transforms related to the wave equation (e.g. [( 2 )]) (see also [( 1 )] and [( 172 )] for related results). However, as the famous counterexample by Boman [( 20 )] shows, the condition of infinite smoothness of the weight function ω alone does not guarantee uniqueness.
Let us move now to the much simpler case of the exponential transform R_{μ}. One of the most important (albeit simple) properties of the Xray transform is the Fourierslice (also called projectionslice) theorem (see Chapter 2 and (A.3)). There is an analog of this property for the exponential transform [( 121 ), ( 124 )], Namely, let f be a compactly supported function on the plane. A straightforward calculation shows that
Let us turn now to the problem of uniqueness of reconstruction, i.e. injectivity of the operator R _{μ}. We will follow here considerations of [( 117 )] and [( 74 )]. Let us assume that function f is either compactly supported or decays as exp(−(μ +∈) x) for a positive ∈. Consider the surface
As soon as injectivity is established, it is natural to look for inversion formulas. By now many inversion formulas have been established for the attenuated Radon transform (see [( 8 ), ( 128 ), ( 129 )] and also [( 26 ), ( 43 ), ( 58 ), ( 91 ), ( 96 ), ( 98 ), ( 125 )] for further discussion and numerical implementation). We will address these below in this section. We would like to start with the much simpler, but still very rich, example of the exponential transform.
A.3.2.1 Exponential transform
The first explicit inversion formula for the exponential Xray transform in the plane was provided by Tretiak and Metz in [( 167 )] (see also its discussion in [( 124 )]). An inversion procedure was also provided in [( 12 )]. Let us introduce the dual exponential Xray transform (or exponential backprojection) ${R}_{\mu}^{\#}$ as follows: applied to a function g(s, ω) it produces a function on the plane according to
Theorem A.5 [( 167 )] Let $f\in {C}_{0}^{\infty}\left({\mathbf{R}}^{2}\right)$. Then
The original proof of this filtered backprojectiontype formula in [( 167 )] and the proof provided in [( 124 )] are interesting and instructive, although rather (p.699) technical. We will now present a different simple approach based on complex analysis [( 89 ), ( 95 )].
Let us go back to the projectionslice formula (A.12). As we have already discussed, it implies that given R _{μ} f one can obtain the values of the Fourier transform of the unknown function f on the surface S _{μ} in C ^{2}. If instead of S _{μ} we dealt with R ^{2}, we would just use the Fourier inversion formula to find f. So, the question is whether one can develop a Fourier inversion formula that uses the data on S _{μ} instead of R ^{2}. The simple idea is to try to use the Cauchy theorem in order to deform the surface of integration from R ^{2} to S _{μ}. Let us look at the standard Fourier inversion formula in R ^{2}:
(p.700) It is contained in the complex line z _{2} = iz_{1} and hence is analytic. We now define the following part of the surface S _{μ}:
The idea of using the Cauchy theorem and Cauchy formula for inverting the exponential Xray transform was initially used in papers [( 44 ), ( 117 ), ( 121 )], and developed in the presented form in [( 89 ), ( 95 )]. This approach leads to a wide variety of inversion formulas [( 95 ), ( 158 )], including cases when the attenuation μ depends on the direction ω [( 95 ), ( 158 )] and when imperfect detector collimation is taken into account [( 89 )]. Invertibility of the exponential Xray transform was also studied analytically and numerically by different methods in [( 63 )−( 65 )]. One might also be interested in related inversion formulas available for the cases of halfview [( 127 )], threedimensional limited data acquisition [( 97 )], and imperfect collimation [( 89 )].
Another interesting relation with complex analysis is the following. The possibility of pasting disks to totally real surfaces in R ^{2} is a popular topic of complex analysis (e.g. [( 55 )]). In particular, the question of rigidity of such a pasted disk has been studied. Rigidity means absence of parametric deformation of the pasted analytic disks. The natural way the disk (A. 14) pasted to the surface S _{φ} appears, suggests that one might expect its rigidity. This happens to be true [( 112 )]. Moreover, a general totally real surface in C ^{2} in a neighborhood of the boundary of a pasted analytic disk can be reduced to a canonical form resembling the definition of S _{μ}, which in turn leads to the possibility of counting the number of parameters in deformations of this disk [( 112 )]. This provides a different approach to the known results [( 55 )] concerning this problem.
A.3.2.2 Attenuated transform
We follow in this section the technique introduced in [( 8 ), ( 25 )] with an additional development from [( 91 )]. Let as before f(x) be the source's intensity in a planar domain ω, u(x, φ) the density of particles at x moving in the direction of ν_{φ} = (cosφ, sinφ) and μ(x) the linear attenuation coefficient. When scattering is absent, the transport equation (A.10) can be rewritten as
(p.701) Using the complex coordinate z = x_{1} − ix_{2} one gets
The boundary data u(x, φ)_{x∈∂ω} known. Let us consider the spaces L _{2}(S) and H ^{2} ⊂ L _{2}(S) (the Hardy space) on the circle and the orthogonal projector P : L _{2}(S) → H ^{2}. We also denote by u _{k} = u ^{ k } (x) the Fourier coefficients of u(x, φ). Let R be the right shift in H ^{2}:(Ru)_{k} = u _{k−1} and L = R* be the left shift. Then we easily obtain from (A. 16) that
Let υ = P _{υ} = υ_{e} ⊕ υ_{0}, where υ_{e} = (υ_{0}, υ_{2},…) and υ_{0} = (υ_{1}, υ_{3},…). Then ( A.18 ) can be rewritten as
The boundary values of both vectors υ_{e} and υ_{0} at ∂Ω are known. According to ( A.17 ), only the odd part υ_{0} plays a role in the recovery of f(x).
Let us introduce a generalized ∂¯ operator:
The equation ( A.19 ) is a generalized ∂¯ equation for the function υ If one could prove uniqueness of its solution under given boundary conditions, this would immediately imply uniqueness of recovery of f, and hence uniqueness for the attenuated Radon transform. This would follow, for instance, from a Cauchytype formula for such “analytic” vectorvalued functions (if A = L ^{2}, such functions are called Aanalytic in [( 8 ),( 25 )]). This is exactly what will be outlined in the next subsection.
(p.702) A.3.2.3 Cauchytype formulas and inversion of Radon transform
Let ω be a planar domain and u be an H ^{2}valued function. Then analogs of the standard considerations lead to a Cauchy–Greentype formula:
Theorem A.6 If u ∈ C ^{ l }(ω¯,H^{2}) and z _{0} ∈ ω then
Corollary A.7 If Du = 0 (i.e. u is Aanalytic), then

(1) Find the boundary values for u _{e} and u _{0}.

(2) Find u _{e} and u _{0} using ( A.21 ).

(3) Find f from υ_{0} using ( A.17 ).
Remark 1. As one can read in [( 8 ), ( 25 )], the resolvent (z − z _{0} − L ^{2} (z¯ − z¯_{0})) involved in the above formulas is formally not defined, since the point λ = (z − z _{0}) / (z¯ − z¯_{0}) always belongs to the spectrum of the operator L ^{2}. It can, however, be extended from the exterior of the disk λ ≤ 1 as a strong limit and hence defined as an unbounded operator.
Remark 2. It is interesting to recognize the standard parts of the Radon inversion formula (Hilbert transform, differentiation, and back projection) in the procedure above. The continuation of the resolvent up to the spectrum gives the Hilbert transform, the integration in the Cauchy formula ( A.21 ) corresponds to the back projection, and the differentiation is performed in ( A.17 ).
A.3.2.4 Inversion of the attenuated transform
We have dealt so far with the case of zero attenuation. Let us outline now how one can resolve the situation with the attenuated transform, i.e. when the transport equation has the form
We start with the case of the exponential transform, which is simpler than the one of the general attenuated transform and also leads to interesting analysis.
A.3.3.1 Exponential transform
Probably the first appearance of the range conditions for this transform was in the papers [( 12 ), ( 167 )] devoted to its inversion. We consider a compactly supported smooth function f(x) in the plane and its exponential Xray (Radon) transform with attenuation μ:
(p.704) We consider the Fourier transform ∶ι _{l}(σ,μ) of g_{l} (s,μ) with respect to s. It was observed in [( 12 ), ( 167 )] that the function
Theorem A.8 A function g(s, ω) can be represented as R_{μ} f for some $f\in {C}_{0}^{\infty}\left(R\right)$ if and only if

(i) $f\in {C}_{0}^{\infty}(R\times {S}^{1})$

(ii) function (σ +μ)^{ι}ĝ(σ, μ) is even with respect to σ for any ι ∈ Z.
An interesting corollary of the range condition provided in this theorem is that the function ĝι(σ, μ) has a zero of order ι at the point μ × sgn ι. Existence of these zeros, however, is obviously not equivalent to the whole set of conditions. One can check directly that ( A.30 ) in the case of constant attenuation leads to the abovementioned root conditions. This shows in particular that the range conditions ( A.30 ) for the attenuated transform that will be discussed later are incomplete [( 93 )].
The range conditions of Theorem A.8 do not have the usual momentum form as they do in the case of the Radon transform. A momentumtype set of conditions was also found in [( 93 )] and [( 94 )].
Theorem A.9 A function g(s, ω) can be represented as R_{μ}f for some $f\in {C}_{0}^{\infty}\left({R}^{2}\right)$ if and only if

(i) $f\in {C}_{0}^{\infty}(R\times {S}^{1})$

(ii) the following identity is satisfied for any odd natural n:
(A.24)Here i is the imaginary unit, ω = (cosφ, sinφ), and ο denotes composition of differential operators.$$\mathrm{\sum}_{k=0}^{n}\left(\begin{array}{c}n\\ k\end{array}\right)\frac{d}{d\phi}o(\frac{d}{d\phi}i)o\dots o(\frac{d}{d\phi}(k1\left)i\right){\int}_{\infty}^{\infty}\left(\mathrm{\mu s}{)}^{nk)}g\right(s,\omega )\mathrm{ds}=0.$$
The condition ( A.24 ) is not very intuitive. One way we can try to understand it is the following. Let us assume that g is representable as R _{μ}f and plug R _{μ} f instead of g into ( A.24 ). In the case of the standard Radon transform, this would lead to a condition that is obviously satisfied. In other words, for the Radon transform, necessity of the momentum conditions (as soon as they are formulated) is simple to observe, and only sufficiency requires a proof. However, even necessity of ( A.24 ) is not obvious. A direct calculation shows that necessity of this condition is equivalent to the following series of identities for the function sinφ.
(p.705) For any odd natural n
A better understanding of the meaning of the range conditions for the exponential Radon transform can be achieved using complex analysis. This was done in the series of papers [( 3 ), ( 4 ), ( 131 ), ( 132 )]. It was shown that the range description problem essentially reduces to BernsteinHartogs’ type of theorems on extension of separately analytic functions (see Chapter 5 and further details later on in this chapter).
Consider again formula ( A. 12 ):
Returning to the projectionslice formula ( A.12 ), we notice that although we initially used the real values of σ only, we can in fact use all complex values of σ provided that the function f is compactly supported. This leads to the determination of values of f̃ on a threedimensional variety in C ^{2}. Let us denote this variety by ∑_{μ}:
Theorem A. 10 Let F be a function defined in the exterior of the disk of a positive radius in the plane. If this function is entire along each tangent line to the disk, it is extendable to an entire function of two variables.
This theorem was proven (although not formulated explicitly in the present form) in [( 3 )]. The proof uses an analytic trick, whose role has not yet been well understood. For instance, the question remains whether one can replace tangent lines to a circle by tangent lines to a more general convex algebraic curve. This is certainly a separate analyticity theorem of Bernstein's kind [( 5 ), ( 18 )], since in a neighborhood of any point in the exterior of the disk the pairs of tangent lines provide a local coordinate system, for which we have separate analyticity. A proof that uses purely complex analytic tools was provided in [( 132 )]. It was shown in [( 136 )] and [( 131 )] that the strange moment conditions ( A.24 ) can be understood as the infinitesimal form of ( A.28 ) at the boundary of the disk.
The last theorem leads to the following range description [( 3 )]:
Theorem A.11 Let g(s, ω) be a function on R × S ^{1}, and μ be a positive real. Then g = R _{μ} f for some function $f\in {C}_{0}^{\infty}\left({R}^{2}\right)$ if and only if $g\in {C}_{0}^{\infty}(R\times {S}^{1}),$ and the Fourier transform ĝ(ξ,φ) with respect to the s variable of g(s, ω(φ)) satisfies the following condition for all real σ:
Here ω(φ) = (cosφ, sinφ).
(p.707) The condition ( A.29 ) is a simple restatement of ( A.28 ).
A generalization of this result to the multidimensional case is provided in [( 4 )]. A nice discussion can be found in [( 131 )] and [( 132 )].
Let us mention briefly some applications of these range conditions. It was shown in [( 137 )] how the conditions can be effectively used in SPECT in order to detect and correct some data errors arising from hardware imperfection. The paper [( 111 )] uses the range conditions in order to treat numerically incomplete data problems of SPECT. It is suggested in paper [( 32 )] that the range conditions might be used in some problems of radiation treatment planning (see also [( 31 )] and [( 88 )]). We will not go into details of these results.
A.3.3.2 Attenuated transform
The question about the range descriptions for the attenuated Radon transform turns out to be a tricky one (see also Section A.2 ). In fact, one could think that there is no chance of getting any explicit range conditions for the attenuated Xray transform. Despite this natural belief, Natterer [( 124 )] found the following infinite set of such range conditions. Let us introduce some notation first. We denote by H the Hilbert transform on the line
Theorem A. 12 [( 124 )] Let f and μ belong to the Schwartz space S(R ^{2}). Then, for k > m ≥ 0 integers, we have
As has already been explained, the set of conditions provided in this theorem is incomplete. Nevertheless, it is of practical importance, since it can be used for the simultaneous recovery of the sources f(x) and attenuation μ (x) (see [( 121 )−( 123 )] and further discussion below).
The question had remained on how to write down the missing conditions, until the very recent publication [( 130 )], where a complete set of conditions was obtained.
It can be shown [( 91 )] that the approach of [( 8 )] provides an explanation and natural derivation of ( A.30 ). In order to explain this, let us start with the case of zero attenuation (i.e. standard Radon transform) and look at it in terms of .Aanalytic functions (see Section ( A.3.2.2 )). In the notation of that section, the role of the Radon data is played by a function u on ∂ω with values in the Hardy (p.708) space. It belongs to the range iff it is 4analytically extendable to ω. Due to an analog of Sokhotsky's formula for such functions, it is easy to show that this happens iff
How can one now include attenuation in this scheme? Using an appropriate “correction” function h as in Section ( A.3.2.2 ) and then applying the above conditions to hu one obtains the range conditions. In particular, one can obtain ( A.30 ).
(p.709) Theorem A.14 [( 91 )] Orthogonality to the suitably understood functions H ^{ k } h (i.e. analogs of orthogonality to z ^{ k }) gives Natterer's range conditions ( A.30 ).
We would like to mention again that a full set of range conditions was obtained in [( 130 )].
As we have discussed in the beginning of this section, simultaneous recovery of the source's density f(x) and of the attenuation μ (x) is an important applied issue. At the first glance this problem might look hopeless. Indeed, in the plane case the data g = T_{μ}f is a single function of two variables, while we are trying to recover two functions f(x) and μ (x) of two variables. However, this counting of functions and variables would be persuasive only if the operator T_{μ} were close to an invertible one (for instance, if it were a Fredholm operator). In fact, the study of the range conditions shows that the operator T_{μ} has an infinite dimensional cokernel. This means that when μ changes, the range could in principle rotate in such a way that for two different values of the attenuation the two ranges would have zero (or a “very small”) intersection. If this were true, then both f and μ would be recoverable or “almost recoverable.” For instance, some additional a priori information could help. We will describe now some results that show that something like this does happen, although the understanding of the situation is still very superficial.
Let us start with the simplest case of the exponential Xray transform. The first attempt to recover both f and μ from R_{μ}f was made in [( 75 )]. This problem was resolved in [( 162 )] (see also [( 163 )]). Consider the range conditions in the form of evenness of the function
(p.710) Theorem A.15 [( 162 ), ( 163 )] Simultaneous recovery of a compactly supported function f and constant attenuation μ from the values of R _{μ} f is possible if and only if f is not radial.
In fact, in the nonradial case one can provide an algorithm of finding the attenuation μ (see [( 162 )] and [( 163 )]).
An alternative approach to recovery of attenuation was described in [( 7 )].
Recovery of a variable attenuation is definitely much more difficult. As in the exponential case, the range theorems are used to this end. Theorem A.12 has been used under some additional restrictions for the sources distribution in order to accomplish this task [( 121 )−( 123 )]. Papers [( 8 ), ( 128 )−( 130 )] contain some additional indications on how their techniques and results might be employed for that purpose.
Electrical impedance tomography (EIT) is a very promising, inexpensive, and technologically simple (at least in comparison with Xray, SPECT, PET, and MRI scans) method of medical diagnostics and of industrial nondestructive testing (see, for instance, [( 9 ), ( 10 ), ( 16 ), ( 17 ), ( 23 ), ( 28 ), ( 29 ), ( 47 ), ( 54 ), ( 77 ), ( 87 ), ( 104 ), ( 150 ), ( 154 ), ( 155 ), ( 164 ), ( 165 )] and references therein). Here is the idea of EIT: one places electrodes around a body, creates some known currents through the electrodes, and then measures the corresponding boundary voltage drops. All of the data is collected on the boundary. After many such measurements one wants to recover the distribution of the electric conductivity inside the body. The information about the electric conductivity is very important for medical diagnostics; it is also vital for some electrical procedures, such as defibrillation; it could also provide a cheap and reliable nondestructive evaluation technology. There are several strong groups in industry and academe that work on the practical and theoretical aspects of this problem using various techniques. We will consider here only one approach, which is related to the integral geometry problems discussed in the book. The reader can refer to the papers quoted above that address other methods that do not use the Radon transform.
Let us describe first the setup of the inverse conductivity problem, which is the mathematical formulation of EIT. Let U ⊂ R^{ n } be a domain with sufficiently smooth boundary γ. An unknown function β(x) (the conductivity) must be recovered from the following data. Given a function ψ on γ (the current) one solves the Neumann boundary value problem
(p.711) In other words, the socalled DirichlettoNeumann operator λ_{β} : φ → ψ is known. One needs to recover the conductivity β from this data. The problems is obviously nonlinear.
The inverse conductivity problem is a hard nut to crack both analytically and numerically. The main questions are about uniqueness of determination of the conductivity, stability of the reconstruction, and inversion algorithms. Although the problem of uniqueness can be considered as principally resolved (see [( 78 ), ( 119 ), ( 166 ), ( 170 )], and references therein), the other ones are still under thorough investigation. The general understanding is that the problem is highly unstable, so there is no hope of achieving the quality of reconstruction known to other common tomographic techniques, at least without involving additional information about the image to be reconstructed.
The first practical algorithm of Barber and Brown [( 9 )−( 11 )] deals with the linearized problem, which is still highly unstable. A thorough investigation of this algorithm started by Santosa and Vogelius [( 155 )] lead in the works of Berenstein and Casadio Tarabusi [( 16 ), ( 17 )] to the understanding that the linearized twodimensional problem can be treated by means of hyperbolic geometry (see Section ( 7.3 )). Consider the twodimensional case when U is the unit disk. It is well known (see for instance [( 13 ), ( 69 )]) that the unit disk serves as the Poincare model of the hyperbolic plane H ^{2}. Here are the indications why the hyperbolic geometry might play some role in the (at least linearized) inverse conductivity problem. First of all, the Laplace operator that arises in the linearized problem is invariant with respect to the group of Mobius transformations. Another indication is that if one creates a dipole current through a point on the boundary of U, then the equipotential lines and the current lines form families of geodesies and horocycles in H ^{2}. Following the analysis done in [( 155 )] of the algorithm suggested in [( 9 ), ( 10 )], Berenstein and Casadio Tarabusi [( 16 ), ( 17 )] discovered that in fact for n = 2 the linearized inverse conductivity problem can be reduced to the following integral geometry problem on H ^{2}: the available data enables one to find the function
It seemed at first that the hyperbolic geometry approach would not work in dimensions higher than two, due to lack of hyperbolic invariance of the governing equations. It was shown, however, in [( 46 )] that a combination of euclidean and hyperbolic integral geometries does the trick.
Notes:
(1) The second author was partially supported by NSF grants 9877155 and 0200788.