Almost-everywhere extrapolation from cone-beam data

Abstract

A source applies imaging energy that passes through an object being imaged and is then detected by a detector. A scanning trajectory causes the detector to provide a Cauchy data set of measured cone beam data to a processor. The processor extrapolates by use of John's equation such that missing cone beam data is determined. The measured cone beam data and the determined cone beam data together provide a complete data set for exact image reconstruction.

Description

BACKGROUND OF THE INVENTION

This invention relates to X-ray imaging techniques, and, more particularly, to computerized tomography X-ray imaging using a minimal set of cone beam data measurements.

Computerized tomography (CT) uses a source of imaging energy such as X-ray energy and a detector for detecting imaging energy that has passed through an object of interest, often a patient being imaged for medical purposes. Typically, a single point source is used with an area detector such as an X-ray detector array.

Relative movement between the source, object of interest, and detector is used to collect data for image reconstruction purposes. Usually the source is moved, while the object and detector remain stationary relative to each other.

Typical source trajectories are "1 dimensional manifolds" described by parametric equations of a single variable. The advent of area X-ray detectors permits measurement of a 2 dimensional data set for each source position, since the detector lies in a 2 dimensional plane. Therefore, it is possible to measure a 1+2=3 dimensional data set of line integrals of the object being imaged.

This dimensional count is encouraging for volumetric imaging such as CT because the imaging object lies in three spatial dimensions. More specifically, complete CT data for a family of parallel planes constitutes a 2+1=3 dimensional data set and completely determines the imaging object.

Unfortunately, cone beam data generated by a single point source does not permit such simple reconstruction. A single circular trajectory does not provide a complete data set for Radon reconstruction (i.e., apart from measurement and discretization errors). Therefore, various scan trajectories have been used.

Among scan trajectories are those disclosed in U.S. Pat. No. 5,465,283 issued Nov. 7, 1995 to Tam. That patent discloses, among other information, the use of a scan trajectory of two offset circular scans with a line extending between them.

Generally, the scan trajectories result from relative movement among the source, detector, and object being imaged. Usually, the source is moved in a path about the object to define the scan trajectory. The object being imaged is often a medical patient, but could also be an industrial part being imaged to locate possible defects. When the object being imaged is an industrial part, the scan trajectory may be defined at least partly by movement of the object relative to the source and/or relative to the detector. A scan trajectory might also be defined at least partly by movement of the detector.

Such scan trajectories that provide complete cone beam data may present a number of problems. The path or paths followed by the source or other component being moved may be complex. This requires a complex robotic function. Additionally, the scan time for obtaining a complete data set may be longer than desirable. Especially in the case where the imaged object is a patient, it is desirable to limit the time and dose of X-ray exposure. A complex scan path requires a longer time of exposure than a simple path. Yet, if the dose is reduced to partially compensate for a long exposure time, the signal to noise ratio is reduced. Errors in the data, particularly due to lag and signal to noise degradation with decreasing dose, often are more troublesome with complex scan paths. Further, a complex scan path increases the amount of measured data and this in turn increases computational demands when using the reconstruction process to provide an image.

Although it is possible to provide image reconstruction with a incomplete data set, such a reconstructed image will be quite inexact, if not undetermined, in some regions as a result of not having the missing data.

BRIEF SUMMARY OF THE INVENTION

In an exemplary embodiment of the method of the invention, an object is imaged by applying imaging energy from a source to the object. Imaging energy that has passed through the object is detected by a detector. The object is scanned with the imaging energy such that the detector collects measured image data that constitutes a Cauchy data set for John's equation. Cauchy data sets and John's equations are known mathematic concepts. Briefly, a Cauchy data set allows one to solve a differential equation in a particular region enclosed by a Cauchy surface. John's equation sets consistency conditions for line integrals passing through the object, the line integrals corresponding to imaging energy passing through the object. Extrapolation is performed on the Cauchy data set to solve John's equation and determine missing cone beam data.

An exemplary embodiment of the system of the invention includes a source for applying imaging energy to the object. A detector detects imaging energy that has passed through the object. A positioner scans the object with the imaging energy such that the detector collects measured image data that constitutes a Cauchy data set for John's equation. An extrapolator extrapolates from the Cauchy data set to solve John's equation and determine missing cone beam data. An image supplier supplies an image of the object based on the measured image data and the determined missing cone beam data.

This technique has a number of advantages over the prior art. The scanning time (data acquisition time) and reconstruction time can be reduced. Errors in the data, particularly due to lag and signal to noise degradation with decreasing dose, can be minimized. The amount of measured data is minimized and this reduces the computational demands for image reconstruction process steps.

BRIEF DESCRIPTION OF THE DRAWINGS

FIG. 1 is a simplified illustration of components of an imaging system and showing process steps within a processor;

FIG. 2 is an illustration of some components of the imaging system and showing line integrals measured by the system;

FIG. 3 shows a graph to illustrate mathematic concepts;

FIG. 4 shows the geometry of a particular source trajectory; and

FIG. 5 is a planar view of the source trajectory of FIG. 5 illustrating locations corresponding to extrapolated data.

DETAILED DESCRIPTION OF THE INVENTION

Referring now to FIG. 1, imaging energy is applied by source 10 to object 12 and detected by detector 14. The source 10 can be a cone beam X-ray source, whereas the detector 14 can be an area detector such as a two dimensional array detector having an array of individual detector elements (separate elements not shown). The object 12 being imaged can be a patient, an industrial part, or other object being imaged using CT.

Cone beam energy that has passed through the object 12 is converted to corresponding electrical signals and sent to a data acquisition unit 16 that registers the electrical signals. Unit 16 in turn sends cone beam data to a processor 18, which can be a computer programmed to perform various process steps.

Various process steps performed in processor 18 are shown in simplified form in FIG. 1. At block 20, the measured image data provided is not a complete set for image reconstruction. However, the measured image data is a Cauchy data set for John's equation. A more detailed explanation of the significance of Cauchy data sets and John's equation follows below. At this stage, the Cauchy data set of measured data allows one to calculate missing cone beam data by extrapolation performed by the processor 18 at block 22. The measured data of block 20 and the missing cone beam data determined at block 22 together constitute a complete data set for exact reconstruction (i.e., meaning subject only to measurement errors and discretization errors). This complete data set of cone beam data, some measured and some extrapolated, is supplied to block 24 where processor 18 processes it in the usual fashion for cone beam data. For example, Radon derivative data is calculated followed by calculation of Radon data. Inverse Radon transformation may then be performed to provide reconstructed image data at block 26. The processing of the complete cone beam data may be, for example, as disclosed in the above listed Tam patent, which patent is hereby incorporated by reference. Regardless of whether that or another image reconstruction process from cone beam data is used, the block 26 operation of processor 18 supplies an image of the object (i.e., data corresponding to an image) that may be displayed by display 28, stored in a memory 30, or both.

The scanning trajectory will be defined by relative movements of the source 10, object 12, and detector 14 respectively controlled by robotic positioners 32, 34, and 36. Each of the robotic positioners 32, 34, and 36 causes relative scanning movement among the source 10, object 12, and detector 14. The positioners 32, 34, and 36 are in turn controlled by processor 18 to realize a scan trajectory satisfying criteria discussed below.

Although two or all three of the robotic positioners 32, 34, and 36 might be part of the imaging system of FIG. 1, only one of the three would usually be in a given system. For example, the positioner 32 could be used to determine a scan trajectory by moving source 10 relative to object 12 and detector 14, both of which remain stationary. In that example, there would not be positioners 34 and 36.

FIG. 2 will be used to explain basic concepts used by the present invention before presenting more complex mathematical explanations. A source 10 is shown applying imaging energy beams 38 to pass through object 12 and strike the detector 14 (shown in side view). Four such energy beams 38 are shown for ease of illustration, but as readily appreciated the actual output of source 38 is a cone beam which may be considered as individual rays corresponding to the resolution of the area detector 14. In the simplified illustration, the detector 14 detects data corresponding to the four rays 38. Specifically, the data detected by detector 14 is dependent on line integrals 40 representative of the extent to which rays 38 are attenuated by passing through object 12.

One could measure a complete data set for exact reconstruction (i.e., providing a three-dimensional image) by measuring a sufficient set of line integrals 40 as source 10 is scanned in a trajectory. However, scanning in such a trajectory is complex and subject to other disadvantages.

Instead of measuring a complete data set by a complex scanning trajectory, the present technique is to obtain a data set corresponding to line integrals such as 40 where the data set is at least a minimal set of cone beam measurements such that missing cone beam data can be determined. The line integrals such as 40 must satisfy certain criteria that can be used to solve for missing cone beam data. Given a Cauchy data set for John's equation as explained in more detail below, one can solve John's equation and determine the missing cone beam data. A minimal amount of data for image reconstruction may thus be collected. This provides minimal scan and reconstruction times. Also, this reduces the X-ray dose to a patient when the object being imaged is a patient.

The space of lines in three real dimensions has greater dimension than the space in which the object being imaged sits. The present techniques are designed to take advantage of the overdeterminedness of the map taking an object to its line integrals by extrapolating the X-ray transform of the imaging object from a minimal cone beam data set of the object. This permits exact reconstruction (up to measurement and discretization errors) of the object being imaged. Furthermore, the reconstruction methods used are not predicated upon the specific extrapolation technique described below. The present technique works with a minimal cone beam data set (i.e., no more data than a Cauchy data set), but additional data can be measured to smooth noisy data in some instances such as imaging industrial parts.

Detailed Mathematics

The technique described below is designed to extrapolate almost everywhere the X-ray transform of an object from a minimal set of cone-beam measurements. The detailed mathematics is discussed below in the order of the following outline:

a. definition of the X-ray transform and its relationship to cone-beam data b. comparison of the "parameter counts" of possible density functions of the imaging object to that of the X-ray transform c. derivation of consistency conditions upon the X-ray transform (John's equation) d. derivation of the fundamental solution to John's equation e. extrapolation of data, solution to the Cauchy problem via the fundamental solution to John's equation

Points a, b, and c are simply the backdrop against which point d was developed. The invention uses the fundamental solution to John's equation to solve almost everywhere the Cauchy problem for the "rest" of the X-ray transform, meaning the missing cone beam data.

Definition of the X-ray Transform

After taking -ln (natural logarithm), the data we measure approximates line integrals of the imaging object. There are many ways to parametrize the space of lines in R.sup.3. For our purposes, we first define the X-ray transform in terms of .eta., .xi., .epsilon. R.sup.3, where .eta..noteq..xi.. We define the X-ray transform of a function, f, to be the line integral of f along the line through .xi. and .eta. and scaled as follows: ##EQU1##

After taking -ln, the data that we measure is a normalized X-ray transform, where the direction vector is normalized to have length one. ##EQU2## Parameter Count Comparison

Notice that the X-ray transform can be specified almost everywhere with only four parameters. All lines in R.sup.3 except for horizontal lines can be specified by points .xi. such that .xi..sub.3 .ident.0 and vectors .eta. such that .eta..sub.3 .ident.1. See FIG. 3. X f(.xi.;.eta.) integrates f along the line passing through .xi. and .eta.. Note that almost all lines except for horizontal lines can be parametrized by .xi. lying in the plane defined by (z=0) and .eta. in the plane {z=1}.

This means that the X-ray transform is essentially a function of four variables, .eta..sub.1, .eta..sub.2, .xi..sub.1, and .xi..sub.2. ##EQU3##

For sake of simplicity in deriving John's equation, assume that the function f .epsilon. C.sub.n.sup.2 (R.sup.3). Then ##EQU4##

The function f is dependent upon only three variables, whereas the function X f is written as a function of four variables. The mapping X is therefore overdetermined. This overdetermination is expressed in John's equation, a second order partial differential equation (PDE) satisfied by X f.

Derivation of John's Equation

For f .epsilon. C.sub.0.sup.2 (R.sub.3) we may differentiate X f with respect to .eta..sub.1, .eta..sub.2, .xi..sub.1, and .xi..sub.2. ##EQU5##

John's equation is simply ##EQU6## Note that in this particular case with .epsilon.=0 and .eta..sub.3 =1 only 0.8 is relevant. Generally, however, the same derivation results in a system of three PDEs of the form 0.8. This applies to cone-beam data in that ##EQU7##

In the next section we will derive the fundamental solution to 0.8 and in the subsequent section use this fundamental solution for X .sub.N f given only measurement of X.sub.N f on a .vertline..xi.-n.vertline. noncharacteristic surface of 0.8.

Derivation of the Fundamental Solution to John's Equation

Since John's equation is a constant coefficient PDE, we compute its fundamental solution using Fourier transforms. The fundamental solution satisfies ##EQU8## We take as dual variables x and y, where x, y, .xi., .eta. .epsilon. R2

.xi..about.x and .eta..about.y

The inverse transform of 0.10 is ##EQU9##

To solve for u=u (.xi.; .eta.) we must solve 0.11 for u and take its Fourier transform. Although this is conceptually clear to those trained in Fourier analysis, the details are carried out for clarity below. Before beginning, define

.crclbar..sub.v =angle between vector v .epsilon. R.sup.2 and e.sub.1 =(1,0) ##EQU10##

Using this notation we compute the Fourier transform of u in polar coordinates ##EQU11##

In order to evaluate 0.13 we make more changes of variables

.sigma.=.vertline..xi..vertline. cos (.crclbar..sub..alpha. -.crclbar..sub..xi.) d.sigma.=-.vertline..xi..vertline. sin (.crclbar..sub..alpha. -.crclbar..sub..xi.)d.crclbar..sub..alpha. .crclbar..sub.a =.crclbar..sub..xi. -cos.sup.-1 (.sigma./.vertline..xi..vertline.) .sigma. .epsilon. .vertline..xi..vertline.(-cos (.crclbar..sub..xi.), cos (.crclbar..sub..xi.)) (0.14) and .tau.=.vertline..eta..vertline. cos (.crclbar..sub..beta. -.crclbar..sub..eta.) d.tau.=-.vertline..eta..vertline. sin (.crclbar..sub..beta. -.crclbar..sub..eta.)d.crclbar..sub..beta. .crclbar..sub..beta. =.crclbar..sub..eta. -cos.sup.-1 (.tau./.vertline./.eta..vertline.) .tau. .epsilon. .vertline..eta..vertline.(-cos (.crclbar..sub..eta.), cos (.crclbar..sub..eta.)) (0.15)

Combining 0.13, 0.14, and 0.15, we get: ##EQU12## Note first that the denominator is multiple valued, and secondly that u is inversely proportional to the cross products of the vectors .xi. and .eta.. We choose the sign of the denominator so that ##EQU13##

To verify that u is in fact a fundamental solution of John's equation we must verify that u satisfies John's equation away from the origin and that ##EQU14## These facts can be easily verified, but need not be included here. Extrapolation of Data--Solution to the Cauchy Problem

Given a fundamental solution to a PDE, solution of the Cauchy problem is well documented. In this section, we consider a special, but practical case shown in FIG. 4 satisfying the following criterion:

1. the imaging object has finite size and lies between the planes z=0 and z=1 2. the source moves along the circular path or, more specifically on a path .eta..sub.3 =1, .eta..sup.2/1 +.eta..sup.2/2 =1 3. the cone beam envelopes the imaging object 4. the detector is large enough to detect all radiation which has passed through the object (it may lie in the plane .xi..sub.3 =0) For simplicity of notation, write for .xi., .eta. .epsilon. R.sup.2 ##EQU15## and recall that we measure X.sub.N and can easily compute .vertline..xi.-.eta..vertline.. Therefore, it is sufficient to extrapolate d(.xi.;.eta.). Under these conditions we can extrapolate d in .OMEGA. from the data measured on the boundary d.OMEGA. if d.OMEGA. is noncharacteristic. (See equation 1.3) First define: ##EQU16##

We can rewrite John's equation in this notation as

1/2.gradient..multidot.D(.alpha.,.beta.)=d.sub..beta..sbsb.2.sub.,.alpha..sbsb.1 -d.sub..beta..sbsb.1.sub.,.alpha..sbsb.2 =0

For (.eta., .xi.) .epsilon. .OMEGA. we can compute d(.xi.;.eta.) as follows: ##EQU17##

This can be simplified further, by noting: ##EQU18## Plugging these identities into the second integral in 6.1 we get ##EQU19## Plugging 6.3 and the definition of u into 6.1 we get the solution to the Cauchy problem. For (.xi., .eta.) .epsilon. .OMEGA., d(.xi., .eta.) is a weighted integral of d and its derivatimes on .differential..OMEGA.. ##EQU20##

FIG. 4 shows the geometry of the particular source trajectory for which John's equation is derived and data can be extrapolated.

FIG. 5 shows source positions and corresponding extrapolated data corresponding to FIG. 4. This is a view of the plane in which the source trajectory lies. The source moves along a circular path. Extrapolated data is data that could have also been measured from source positions within the circular source trajectory. Note that data extrapolation is not predicated upon a circular source trajectory.

Suitable Cauchy surfaces should be identified, not only for the well-posed boundary value problem, but also for the practical constraints that come with data generated using a point source and measured by a flat panel area detector. The symbol for John's equation is:

.LAMBDA.(.xi.,.eta.)=.xi..sub.1 .eta..sub.2 -.xi..sub.2 .eta..sub.1 Equation (0.21)

If the function .phi.=.phi.(.epsilon., .eta.) define the surface .differential..OMEGA., then the criteria for .differential..OMEGA. to be noncharacteristic is that

.LAMBDA.(.gradient..phi.).noteq.0 Equation (0.22) When .phi.(.epsilon., .eta.)=.epsilon..sub.1 .eta..sub.2 -.epsilon..sub.2 .eta..sub.1 -1 then .LAMBDA.(.gradient..phi.).ident.1 Equation (0.23) If this formulation for an can be used in equation 0.20 (despite the fact that .OMEGA. is not a domain in .sup.4, the solution for d(0,0) is: ##EQU21##

This solution to the Dirichlet problem does not require any derivatives of the measured data. The solution to the more general Cauchy problem, however, does require differentiation of the data.

Consistency conditions upon cone beam data can be derived for any source trajectory and can be solved for any trajectory generating Cauchy data on a non-characteristic surface. To date, the simplest geometry found for data extrapolation appears to be that shown above. However, the technique extends to more general geometries. In the case of an unbounded object (in the CT sense) which is not enveloped by the cone beam, additional measures will be required. It may well be sufficient to augment the missing line integrals corresponding to source positions on the source trajectory with data satisfying integral conditions upon Cauchy data. This data will not necessarily be correct for the given imaging object, but should provide a physical data set which can be used by the extrapolation procedure described above. As the errors in resulting reconstructions can lie as far as the field of view as do the incorrect Cauchy data/line integrals, such errors are acceptable in some situations.

Claims

1. A method of imaging an object comprising the steps of:

applying imaging energy from a source to the object;

detecting imaging energy that has passed through the object by use of a detector;

scanning the object with the imaging energy such that the detector collects measured image data that constitutes a Cauchy data set for John's equation;

extrapolating from the Cauchy data set to solve John's equation and determine missing cone beam data; and

supplying an image of the object based on the measured image data and the determined missing cone beam data.

2. The method of claim 1 wherein the image of the object is an exact reconstruction of the object up to measurement and discretization errors.

3. The method of claim 2 wherein the measured image data and the determined missing cone beam data together constitute a complete data set for image reconstruction purposes.

4. The method of claim 3 wherein the detector is an area detector and the source is a point source.

5. The method of claim 4 wherein the scanning step is accomplished by moving the source in a scanning trajectory relative to the object.

6. The method of claim 4 wherein the scanning step is accomplished by moving the object in a scanning trajectory relative to the source.

7. The method of claim 4 wherein the scanning step uses a circular scanning trajectory and the missing cone beam data determined by the extrapolating step includes data corresponding to points interior to the circular scanning trajectory.

8. The method of claim 4 wherein the scanning step uses a scanning trajectory with some portions lying in a plane and the missing cone beam data determined by the extrapolating step includes data corresponding to points in the plane.

9. The method of claim 1 wherein the measured image data and the determined missing cone beam data together constitute a complete data set for image reconstruction purposes, the detector is an area detector, and the source is a point source.

10. The method of claim 1 wherein the extrapolation step determines missing data from within a region based on data measured at the boundary of the region.

11. A system for imaging an object comprising:

a source for applying imaging energy to the object;

a detector for detecting imaging energy that has passed through the object;

a positioner scanning the object with the imaging energy such that the detector collects measured image data that constitutes a Cauchy data set for John's equation;

an extrapolator extrapolating from the Cauchy data set to solve John's equation and determine missing cone beam data; and

an image supplier supplying an image of the object based on the measured image data and the determined missing cone beam data.

12. The system of claim 11 wherein the image of the object is an exact reconstruction of the object up to measurement and discretization errors.

13. The system of claim 12 wherein the measured image data and the determined missing cone beam data together constitute a complete data set for image reconstruction purposes.

14. The system of claim 13 wherein the detector is an area detector and the source is a point source.

15. The system of claim 14 wherein the positioner moves the source in a scanning trajectory relative to the object.

16. The system of claim 14 wherein the positioner moves the object in a scanning trajectory relative to the source.

17. The system of claim 14 wherein the positioner provides a circular scanning trajectory and the extrapolator determines missing cone beam data including data corresponding to points interior to the circular scanning trajectory.

18. The system of claim 14 wherein the positioner provides a scanning trajectory with some portions lying in a plane and the missing cone beam data determined by the extrapolator includes data corresponding to points in the plane.

19. The system of claim 11 wherein the measured image data and the determined missing cone beam data together constitute a complete data set for image reconstruction purposes, the detector is an area detector, and the source is a point source.

20. The system of claim 11 wherein the extrapolator determines missing data from within a region based on data measured at the boundary of the region.