WO2007115200A2

WO2007115200A2 - Method and system for adaptive and robust detection

Info

Publication number: WO2007115200A2
Application number: PCT/US2007/065680
Authority: WO
Inventors: Yao Xie; Bin Guo; Jian Li; Geng Ku; Lihong Wang
Original assignee: University Of Florida Research Foundation, Inc.
Priority date: 2006-03-31
Filing date: 2007-03-30
Publication date: 2007-10-11
Also published as: WO2007115200A3

Abstract

A method of and system of adaptive and robust detection for use in thermoacoustic tomography is provided. The method includes estimating (410) a bipolar acoustic pulse generated by a focal point from multiple thermoacoustic signals using Robust Capon Beamforming, searching (430) for two peaks of the bipolar acoustic pulse within an interval around a calculated arrival time, wherein the interval is determined from an estimated variance of arrival time differences; and calculating (450) a response intensity at the focal point by calculating a peak-to-peak difference of the two peaks.

Description

METHOD AND SYSTEM FOR ADAPTIVE AND ROBUST DETECTION

Inventors

Yao Xie

Bin Guo

Jian Li

Geng Ku

Lihong V. Wang

FIELD OF THE INVENTION

[0001] The present invention relates to the field of imaging, and more particularly, to image reconstruction.

BACKGROUND

[0002] Thermoacoustic tomography (TAT) applies to a wide span of biomedical applications. The physical basis of TAT lies in the contrast of the radiation absorption rate among different biological tissues. Due to the thermoacoustic effect, when a short electromagnetic pulse, such as microwave or laser, is absorbed by a tissue, the heating results in expansion of tissue that generates acoustic signals. In TAT, an image of the tissue absorption properties can be reconstructed from the recorded thermoacoustic signals. The image reveals the physiological and pathological state of the tissue. TAT has been used in many applications including breast cancer detection. [0003] TAT provides both enhanced imaging resolution and good spatial contrast properties. In current practice, microwave imaging offers excellent contrast between cancerous and normal breast tissue, but the relatively long wavelengths of microwaves provide poor spatial resolution. In contrast, conventional ultrasound imaging has very fine spatial resolution, but poor soft tissue contrast. Accordingly, microwave-induced TAT takes advantage of both the high contrast of biological tissue in electromagnetic frequency band (microwave) and the millimeter range spatial resolution due to the acoustic signals (ultrasound) used in image reconstruction.

[0004] Various image reconstruction algorithms have been developed for TAT. However, developing accurate and robust image reconstruction methods is one of the key challenges encountered in TAT. Approximate reconstruction algorithms, such as the time-domain delay-and-sum beamforming method and the optimal statistical approach, have been proposed. The delay-And-sum (DAS), weighted or unweighted, is a data- independent approach widely used in medical imaging. Image reconstruction using DAS requires minimal prior information on the tissue for image reconstruction. Consequently, DAS can be a fast and straightforward technique for processing wideband acoustic signals. Although not based on an exact solution, DAS algorithms can provide image qualities comparable to those of exact reconstruction algorithms. However, these data- independent methods tend to suffer from poor resolution, high sidelobe level problems, and poor interference rejection capabilities. These image reconstruction algorithms assume that the surrounding tissue is acoustically homogeneous. This approximation may be inadequate in many medical imaging applications. The heterogeneous acoustic properties of biological tissues cause amplitude and phase distortions in the recorded acoustic signals, which can result in significant degradations in image quality. [0005] In ultra-sound tomography (UT), wave-front distortion due to heterogeneity of biological tissue has been well studied with regard to analysis of ultrasonic pulse amplitude distortion and arrival time distortion. Various wave-front correction methods have been proposed to compensate for the distortions. However, the current methods for UT are not highly effective at correcting severe amplitude distortions, and usually involve complicated procedures. Consequently, the methods do not apply well for TAT. Distortion issues in TAT are different from that in UT. In a breast UT, for example, amplitude distortion caused by refraction can be more problematic than phase distortion induced by acoustic speed variation. In TAT, however, even for the biological tissue, such as the breast tissue, with a relatively weak heterogeneity, phase distortion dominates amplitude distortion.

[0006] A need therefore exists for a system capable of mitigating phase and amplitude distortions in microwave-TAT for increasing image resolution and contrast.

SUMMARY

[0007] Data-adaptive approaches, such as the robust Capon beamforming (RCB) method, can have much better resolution and much better interference rejection capability than data-independent algorithms such as DAS. By allowing uncertainties in a receiving array steering vector, adaptive and robust techniques (ART) can be used to mitigate the amplitude and phase distortion problems encountered in TAT. [0008] The invention relates generally to a system and method for adaptive and robust detection of thermoacoustic signals for high contrast and high resolution imaging in TAT. The signal processing principles of operation can be applied to other forms of imaging such as ultrasound imaging and microwave imaging, and are not limited to thermoacoustic imaging. One embodiment of the invention is a method that can include receiving multiple thermoacoustic signals emitted from a focal point in response to a heating of a tissue due to radiation absorption. The source of radiation can be a microwave device, laser, or electromagnetic stimulation device. The thermoacoustic signal can include a pulse component and a noise component, both generated as a result of radiation absorption. A waveform of the pulse component can be estimated by applying RCB to the multiple thermoacoustic signals. ARTs applied within the RCB can mitigate amplitude and phase distortion due to the noise component, which can include multi-path interference and variation in acoustic speed through the tissue. A response intensity can be evaluated from the waveform for reconstructing an image of the thermoacoustic response intensity corresponding to an absorption property of the tissue.

[0009] The method can include measuring an amplitude distortion from the thermoacoustic signals, and allowing an uncertainty in the amplitude distortion of the waveform to accommodate an uncertainty range of the phase distortion. For example, the amplitude distortion may be associated with a refraction property of the tissue. The phase distortion may be associated with an sound speed variation due to the tissue. The method can further include searching for a peak in the uncertainty range to estimate a phase distortion, and comparing a location of the peak to an expected location. A phase distortion correction can be applied to compensate the waveform for the phase distortion by shifting the waveform by the determined phase distortion. [00010] The method also can include enhancing a contrast of the image by using a peak-to-peak difference for the response intensity of the focal point. For example, the waveform can be bipolar, having a positive peak corresponding to a compression pulse and a negative peak corresponding to a rarefaction pulse. The response intensity can be the peak-to-peak difference of the bipolar waveform at a particular peak location. The focal point can be scanned from one or more locations to cover an entire cross- section of the tissue, wherein the scanning can occur at multiple heights. [00011] The response intensity of a focal point can be determined, according to one embodiment of the invention, as follows: detecting a peak; calculating an arrival time of the pulse component based on the peak; aligning one or more thermoacoustic signals based on the arrival time; and compensating for a loss in amplitude due to a propagation decay associated with the tissue. According to one aspect, the searching can be performed within an interval around the arrival time, and the interval can be determined by a difference between a true arrival time and the actual arrival time. According to another aspect, the search range can be based on an estimated variance of arrival time differences. The arrival time can account for homogeneous acoustic properties of the tissue.

[00012] Other embodiments of the invention relate to a method of adaptive and robust detection for use in thermoacoustic tomography. A method according to one embodiment can include estimating an acoustic pulse generated by a focal point from one or more thermoacoustic signals using RCB, searching for a peak in the acoustic pulse for mitigating a phase aberration, and calculating a response intensity from the peak for forming one portion of an image. An array steering vector can be estimated using covariance fitting. An optimal beamformer weight vector can be obtained from the estimated steering vector for estimating a waveform of the pulse. The array steering vector can maximize a power of the pulse subject to a dual constraint optimization that uses the optimal beamformer weight vector.

[00013] A first constraint of the dual optimization can be a function of a covariance of the waveform, a power of the pulse in the waveform, and an array steering vector weighting the pulse. A second constraint of the dual optimization can be an uncertainty set that is least one of a spherical uncertainty and an elliptical uncertainty. The uncertainty set allows the array steering vector to accommodate an uncertainty range during a phase measurement. A peak can be searched within the extended uncertainty range for mitigating phase distortions. The error bound can be a user-provided parameter for adjusting a resolution and a contrast of the image. For example, a smaller error can provide an increased ability of the RCB to suppress interference that is close to the pulse for increasing a resolution of the image. A larger error can provide an increased robustness of the RCB to tolerate distortions due to a small sample size of thermoacoustic signals. The uncertainty set can be decreased for images having low resolution, and increased for images distorted by interference.

[00014] A Lagrange multiplier can be employed within the dual constraint optimization to estimate an amplitude steering vector. An Eigenvalue decomposition can be applied to the covariance for creating a monotonically decreasing function of the Lagrange multiplier. Additionally, Newton's method can be applied to solve for the Lagrange multiplier. Scaling ambiguity in the waveform can be eliminated by forcing a norm of the array steering vector to equal a number of transducers employed in capturing the thermoacoustic signals.

[00015] Waveform estimation of the pulse can include applying a weight vector to a pre-processed signal, wherein the weight vector is determined using the estimated array steering vector in a weight vector expression of a Capon beamformer. The pre- processed signal can be a gain-scaled and time-shifted version of the received acoustic pulse. In practice, uncertainty set allows the array steering vector to extend a range for phase measurement. By accommodating uncertainty in the array steering vector, amplitude distortions and phase distortions can be mitigated.

[00016] Another embodiment of the invention is an adaptive and robust detection method for use in thermoacoustic tomography. The method can include estimating a bipolar acoustic pulse generated by a focal point from one or more thermoacoustic signals using robust Capon beamforming and searching for two peaks of the bipolar acoustic pulse within an interval around a calculated arrival time, and the interval is determined from an estimated variance of arrival time differences. The method also can include calculating a response intensity at the focal point by calculating a peak-to-peak difference of the two peaks.

BRIEF DESCRIPTION OF THE DRAWINGS

[0017] There are shown in the drawings, embodiments which are presently preferred. It is expressly noted, however, that the invention is not limited to the precise arrangements and instrumentalities shown.

[0018] FIG. 1A is a schematic diagram of operative features of a thermoacoustic tomography (TAT) imaging system in accordance with an embodiment of the inventive arrangements;

[0019] FIG. 1 B is a plot of acoustic signals generated using a TAT imaging system in accordance with an embodiment of the inventive arrangements; [0020] FIG. 2 is a schematic diagram of a system for adaptive and robust techniques (ART) in TAT in accordance with an embodiment of the inventive arrangements; [0021] FIG. 3 is a schematic diagram of the processor of FIG. 2 in greater detail in accordance with an embodiment of the inventive arrangements; [0022] FIG. 4 is a method of adaptive and robust detection for high contrast and high resolution imaging in TAT in accordance with an embodiment of the inventive arrangements; and

[0023] FIG. 5 further illustrates a method step of FIG. 4 for estimating an acoustic pulse in accordance with an embodiment of the inventive arrangements.

DETAILED DESCRIPTION

[0024] Detailed embodiments of the present method and system are disclosed herein. However, it is to be understood that the disclosed embodiments are merely exemplary, and that the invention can be embodied in various forms. [0025] The terms "a" or "an," as used herein, are defined as one or more than one. The term "plurality," as used herein, is defined as two or more than two. The term "another," as used herein, is defined as at least a second or more. The terms "including" and/or "having," as used herein, are defined as comprising (i.e., open language). The term "coupled," as used herein, is defined as connected, although not necessarily directly, and not necessarily mechanically. The term "suppressing" can be defined as reducing or removing, either partially or completely. The term "processing" can be defined as number of suitable processors, controllers, units, or the like that carry out a pre-programmed or programmed set of instructions.

[0026] The terms "program," "software application," and the like as used herein, are defined as a sequence of instructions designed for execution on a computer system. A program, computer program, or software application may include a subroutine, a function, a procedure, an object method, an object implementation, an executable application, a source code, an object code, a shared library/dynamic load library and/or other sequence of instructions designed for execution on a computer system. [0027] Referring to FIG. 1A, a TAT imaging system 100 according to one embodiment of the invention is shown. The TAT imaging system 100 can be used to create a tissue image from a reconstruction of thermoacoustic signals. The TAT imaging system 100 can emit a stimulating electromagnetic pulse to a tissue 102. The stimulating electromagnetic pulse can be, for example, a laser or a microwave emitting device that can be targeted at the focal point 104. The biological tissue 102 under testing can absorb the stimulating electromagnetic pulse, which can cause a sudden heat change at the focal point 104. Accordingly, an acoustic pulse can be generated at the focal point 104 due to a thermoacoustic expansion effect in the tissue 102. [0028] Due to the thermoacoustic effect, the heating caused by the pulse results in expansion of the tissue, which in turn generates the acoustic signals. An image of the tissue absorption properties can be reconstructed from captured acoustic signals. The image reveals the physiological and pathological state of the tissue. Accordingly, the TAT imaging system 100 can be used in many applications, including for example breast cancer detection. The cancerous breast tissues typically have a two to five times higher microwave absorption rate than their surrounding normal breast tissues, which has been attributed to an increase in the amount of bound water and sodium within malignant cells.

[0029] An ultrasonic transducer array 106 can record the acoustic signals generated from the thermoacoustic effect. The transducer array 106 may be a real aperture array or a synthetic aperture array formed by rotating a sensor around the tissue 102 and recording the acoustic signals at different locations. In the particular example illustrated in FIG. 1 , the number of transducers in the array is M . In the case of a synthetic aperture array, the number of transducer data acquisition locations can be M . Each of the M transducers can be assumed to be omnidirectional. The recorded acoustic signals can be sampled and digitized. Referring additionally to FIG. 1 B, a typical recorded acoustic signal 120 based on data measured from a breast specimen is shown. The pulse 121 is a backscattered signal from the interface (e.g., a bowl holding the breast tissue). The acoustic signal 120 also contains a reflection pulse 122 that is the acoustic pulse generated from the thermoacoustic effect in the tissue 102. The data model for the sampled and digitized acoustic signal 120 recorded by the m th transducer can be given by the following equation: x_m (n) = s_m (n) + e,,,(n), m = !, ^{■ ■ ■} , M. (1 )

where n is the discrete time index, starting from t₀ after the stimulating electromagnetic pulse 121. The scalar s_m(ή) denotes the signal component, which corresponds to the acoustic pulse 122 generated at a focal point 104, and e_m{n) is the residual term, which includes un-modeled noise and interference. For example, interference can be caused by other sources within the tissue such as multi-path propagation or refraction. [0030] An image of thermoacoustic response intensity /(r) , which is directly related to the absorption property of the tissue 102, can be reconstructed from the recorded data set {χ,,,(ή)} . Herein, the (2-D or 3-D) vector r denotes the focal point location coordinate 106. To form an image, the focal point 104 is scanned at location r to cover an entire cross-section of the tissue 102. For example, the transducers can acquire signals at different heights, wherein for each height, a 2-D cross-sectional image can be reconstructed and a 3-D image can be formed from the 2-D images. The homogeneity of the tissue 102 can determine the distortion in the image. For example, the sound speed in human female breast varies widely from 1430 m/s to 1570 m/s around the commonly assumed speed of 1510 m/s. The heterogeneous acoustic properties of biological tissues cause amplitude and phase distortions in the recorded acoustic signals, which can result in significant degradations in image quality. Accordingly, certain uncertainties can be allowed in ART to deal with amplitude and phase distortions caused by the background heterogeneity.

[0031] In practice, the discrete arrival time of the pulse (i.e., for the m th transducer) can be determined approximately as:

The dependence of the arrival time t_m(r) on r can be omitted hereafter for notational simplicity. Here At is the sampling interval, and the 3-D vector r_m denotes the location of the mth transducer. The sound speed v₀ can be chosen to be the average sound speed of the biological tissue under interrogation. The notation |x| denotes the Euclidean norm of x , and [_jμj indicates rounding to the greatest integer less than y .

The second term in (2) represents the time-of-f light between the focal point and the m th transducer.

[0032] The signal components {s_m(ή)}^_n=] are approximately scaled and shifted versions of a nominal waveform s(t) at the source:

exp(-α||r-r II) r - r ,

where a is the attenuation coefficient in Nepers/m. In TAT, the major frequency components of the acoustic signals take a relatively narrow band, and are usually lower than those in UT. The term a can be approximated as a frequency independent constant. [0033] The acoustic signals can be preprocessed to align all the signals from the focal point r and compensate for the loss in amplitude due to propagation decay. For example, let yJri) denote the signal after preprocessing:

y_m in) = exp (a ||r - r_m |) • ||r - r_m | x_m (n + t_m ); (4)

then the received vector data model can be written as:

y(n) = a_os(n) + e(n), n = -N,~ -,N, (5)

where a₀ is the corresponding steering vector, which is approximately equal to a = [l,---,lf , y(ri) = [y_x(n),-- ^■ ,y_M(n)f , e(n) represents the noise and interference term after preprocessing, and 'T' denotes the transpose. The signal of interest x(m) is the acoustic pulse 122, which has been shifted by a time factor t_m and scaled by a gain term. Here the time interval of interests for the signal y(t) can be defined to be from -N to N . Accordingly, N samples can be taken before and after the approximate arrival time given in Equation (2) for the focal point 106 at r . The value of N can be chosen large enough such that the interval from -N to N covers the expected pulse duration. [0034] In practice, both the amplitude and the phase of the acoustic pulse are generally distorted. The phase of the pulse signal can be considered the arrival time of the signal. A major cause for these distortions is the acoustically heterogeneous background. For example, amplitude distortion is mainly due to the interferences caused by multi-path interference, which is inevitable in the heterogeneous medium. Refraction occurs due to acoustic speed mismatch across the tissue interface. Consequently, acoustic pulses arrive at the transducer via different routes and interfere with each other. In contrast, phase distortion is mainly caused by the non-uniform sound speed. For example, in human female breast tissue the sound speed can vary from 1430 m/s to 1570 m/s. Therefore the actual arrival time can fluctuate around the approximately calculated time given in Equation (2). Moreover, an inaccurate estimate of t₀ (t₀ is aligned with the focal point's signal arrival time) and the transducer calibration error may also contribute to the phase distortion. Amplitude and phase distortion can blur the image, raise the image background noise level, lower the values of the object of interest, and consequently decrease the image contrast. The effects of these distortions can be mitigated by allowing a₀ to belong to an uncertainty set centered at a and by considering the signal arriving within the interval from -N to N . [0035] Referring to FIG. 2, a system 200 for adaptive and robust techniques (ART) in thermoacoustic tomography (TAT) is shown. The system can include a receiver 202 for receiving one or more thermoacoustic signals emitted from a focal point in response to a heating of a tissue due to radiation absorption, wherein a thermoacoustic signal includes a pulse component and a noise component, a processor 204 for estimating a waveform of the pulse component by applying robust Capon beamforming to the multiple thermoacoustic signals for mitigating an amplitude and phase distortion of the waveform due to the noise component, and an image reconstruction unit 206 for evaluating a response intensity from the waveform for reconstructing a portion of an image of thermoacoustic response intensity corresponding to an absorption property of the tissue.

[0036] Referring to FIG. 3, the processor 204 can further include an analysis unit 210 for measuring an amplitude distortion from the thermoacoustic signals, wherein the amplitude distortion is associated with a refraction property of the tissue. The processor 204 can further include a peak detector 214 for searching for a peak in the uncertainty range to estimate a phase distortion and comparing a location of the peak to an expected location, and a phase compensator 216 for providing a phase distortion correction to compensate the waveform for the phase distortion by shifting the waveform by the phase distortion. The phase distortion can be associated with an acoustic speed variation due to said tissue.

[0037] Referring to FIG. 4, a method 400 of adaptive and robust detection for use in thermoacoustic tomography is shown. The method 400 is described herein with reference to FIGS. 1-3, but it should be noted that the method 400 can be practiced in any other suitable system or device. Moreover, the steps of the method 400 are not limited to the particular order in which they are presented in FIG. 4. The method can also have a greater number of steps or a fewer number of steps than those shown in FIG. 4. The ART algorithm illustratively comprises three steps: Step 410, robusst Capon beamforming (RCB) for robust waveform estimation; Step 430, peak-searching for phase aberration mitigation; and, Step 450, intensity calculation for forming the final images. These are steps for adaptive and robust detection for producing high contrast and high resolution imaging in thermoacoustic tomography, according to one embodiment.

[0038] At step 410, an acoustic pulse generated by a focal point from one or more thermoacoustic signals can be estimated using RCB. For example, referring back to FIG. 2, the processor 204 can estimate the waveform of the acoustic pulse generated by the focal point at location r , based on the data model in Equation (5). It may appear that the presence of phase distortion has been neglected by using this data model in the first step. However, by allowing a₀ to be uncertain, some phase distortions can be tolerated which causes minimal performance degradation.

[0039] Referring to FIG. 5, a more detailed approach for implementing step 410 is shown. At step 412, an array steering vector can be estimated using covariance fitting. For example, the processor 204 can employ covariance-fitting-based RCB to first estimate the steering vector a₀. Covariance fitting uses the estimated a₀ to obtain an optimal beamformer weight vector for pulse waveform estimation. By assuming that the true steering vector lies in the vicinity of the nominal steering vector a , the following optimization problem can be considered:

maxσ² subject to R - σ²a₀a₀ > 0,

2 Il - Il (6) σ ,a₀ \\ a_Q - a \\ ≤ ε

where A >: 0 means that the matrix A is positive semi-definite, σ² is the power of the signal of interest, and is the sample covariance matrix. The second constraint in Equation (1) is a spherical uncertainty set, though an elliptical uncertainty set can be used if a tighter constraint is desirable. Notably, the uncertainity set accommodates an uncertainty range during phase measurement. The array steering vector of Equation (6) maximizes a power of the pulse subject to the dual constraint optimization.

[0040] The parameter ε in Equation (1) determines the size of the uncertainty set and is a user parameter. To avoid the trivial solution of a₀ = 0 , the method can require that

I— 2 ε < \\Λ (8)

[0041] The smaller the value of ε is, the higher the resolution and the stronger the ability of RCB to suppress an interference that is close to the signal of interest. The larger the value of ε is, the more robust RCB will be to tolerate distortions and small sample size problems caused by calculating R in Equation (2) from a finite number of data vectors or snapshots. When ε approximates M , RCB will perform like DAS. To attain high resolution and to effectively suppress interference, ε should be made as small as possible. In contrast, the smaller the sample size N or the larger the distortions, the larger ε should be. Since the performance of RCB does not depend critically on the choice of ε (provided that it is set to be a "reasonable value"), such qualitative guidelines are usually sufficient for making an appropriate choice for the value of ε .

[0042] At step 414, a Lagrange multiplier can be employed within a dual constraint optimization. By using the Lagrange multiplier method, the solution to (1) is given by he following equation:

, = _ a — r I I + //R ^Λ TJ ¹- a (9)

where I is the identity matrix, μ > 0 is the corresponding Lagrange multiplier that can be solved from the following equation:

[0043] At step 416, an Eigenvalue decomposition can be applied to the covariance for creating a monotonically decreasing function of the Lagrange multiplier. Consider the

Eigendecomposition on the sample covariance matrix R :

R = uru^r, (11)

where the columns of U are the Eigenvectors of R and the diagonal matrix r consists of the corresponding Eigenvalues γ_λ ≥ γ₂ ≥ --- ≥ γ_M . Let b = U^ra , where b_m denotes its m th element. Then Equation (5) can be rewritten as

[0044] At 418, Newton's method can be employed to solve for the Lagrange multiplier. Note that L(μ) is a monotonically decreasing function of μ , with L(O) > ε according to Equation (8) and lim_//→M L(μ) = 0 < ε , which means that μ can be solved efficiently using, for example, Newton's method.

[0045] After obtaining the value of μ , the estimate So of the actual steering vector a₀ is determined by Equation (4). Notably, there is a scaling ambiguity in Equation (1) by treating both the signal power σ² and the steering vector a₀ as unknowns. The ambiguity exists in the sense that (σ²,a₀) and (σ²/c,c^1/2a₀) (for any constant o 0) yields the same term σ²a_oao . To eliminate this ambiguity, the solution So is scaled to make its norm satisfy the following condition:

[0046] Upon completing steps 410-418 for estimating an array steering vector, an optimal beamformer weight vector can be obtained for estimating a waveform of the pulse. At 420, an optimal beamformer weight vector from the estimated array steering vector can be obtained for estimating a waveform of the pulse. For example, to obtain an estimate for the signal waveform s(ri) , a weight vector can be applied to the preprocessed signals {y(ri)}"₌__N . The weight vector is determined by using the estimated steering vector S₀ in the weight vector expression of the standard Capon beamformer:

[0047] Note that Equation (9) has a diagonal loading form, which allows the sample covariance matrix to be rank-deficient. The beamformer output can be written as:

WH) = WRCBy(W), n = -N, ' ~ , N, (15)

which is the waveform estimate for the acoustic pulse generated at the focal point at location r . Upon completion of steps 410-420, the weight vector has been applied to the pre-processed signal for estimating a pulse waveform. The weight vector was determined using the array estimated steering vector in a weight vector expression of a Capon beamformer. Notably, the estimated steering vector was generated by allowing an uncertainty set centered around the norminal steering vector. The pre-processed signal of Equation (4) was a gain-scaled and time-shifted version of the acoustic pulse. [0048] Referring back to FIG. 4, at step 430, the estimated pulse waveform can be searched for a peak for mitigating a phase aberration. This is the second major step in ART. Recall, the estimated pulse waveform was generated in steps 410-418. Briefly, in a homogeneous background, where phase distortion is absent, the arrival time of the acoustic pulse generated by the focal point at location r can be accurately calculate using Equation (2). However, this is rarely true in heterogeneous biological tissues as a result of the non-uniform sound speed through the heterogeneous tissue. When the heterogeneity is weak, such as in the breast tissue, amplitude distortion caused by multi-path is not severe. Accordingly, the original peak can be assumed to remain as a peak in the waveform estimate for the acoustic pulse generated in step 410. [0049] Referring to FIG. 3, the peak detector 214 searches for a peak in an uncertainty range to estimate a phase distortion. The peak detector 214 can search for the peak within an interval around an estimated arrival time, wherein the interval is determined from an estimated variance of arrival time differences. The processor 204 calculates an arrival time of the pulse component based on the peak; aligns one or more thermoacoustic signals based on the arrival time, and compensates for a loss in amplitude due to a propagation decay associated with the tissue. [0050] The bipolar acoustic pulse has one positive peak and a negative peak. The positive and negative peak values can be determined as follows:

P⁺ (18)

P^' (18)

where the searching range [-Δ,Δ] ε [-JV₉N] is around the calculated arrival time given by Equation (2). Here, Δ is a user parameter. Since the peak searching is independent of the particular waveform estimation methods, we use s(ri) to denote the waveform estimated by either DAS or ART. The search range is determined by the difference between the true arrival time F_n, and the calculated arrival time t_m based on Equation

(2). This arrival time difference has been analyzed for breast tissue by taking into account its relatively weak heterogeneity acoustic property . An expression for this difference can be given by:

where r' is a point within the line connecting the focal point at location r and the m th transducer at location r_OT , and v(r') is the local sound speed. In (3), the higher order terms of [v(r')-v_o]/v_o have been ignored. It is reasonable to assume that v(r') is Gaussian distributed with mean v₀ and variance σ². Consequently the arrival time difference is also Gaussian distributed with zero-mean and variance σ] ∞ σ_v ²/v₀ ². If we choose Δ = σ_s , and the duration of the acoustic pulse is τ , we can find the two peaks of the pulse within the interval (-σ_δ,σ_δ +τ) on the recorded signals with a high probability of 0.6826. A symmetric range [-Δ,Δ] around the estimated arrival time can perform similarly to the asymmetric range [-Δ,Δ + r] . The symmetry of the former is easier to handle in practice.

[0051] There can be a tradeoff in choosing the searching range. The larger the searching range, the higher the probability of finding peaks of the acoustic pulse within the range. However, if the range is chosen too large, the interferences may cause false peaks, and as a consequence, the chance of detecting false peaks is greater. In practice, a search range can be empirically based on an estimated variance of the arrival time difference σ_s ■ For example, the search interval can be determined from an estimated variance of arrival time differences, wherein the variance is generated is step 410 by perturbing an array steering vector for allowing an uncertainty in phase measurement.

[0052] At step 450, a response intensity can be calculated from the identified peak for forming one portion of an image. After estimating the waveform generated by the focal point at location r in step 410, a response intensity can be obtained based on the estimated waveform. Two different types of response intensity measurement approaches are amplitude based and energy based, respectively. The waveform peak values obtained in Step 430 of ART can be used for both approaches. Conventional DAS uses the amplitude based measure for TAT imaging, with the corresponding response intensity given by S(O) , or equivalently:

M i_c = s(θ) = ∑y_mΦ), (21)

where the subscript "C" stands for "Conventional." The energy-based measure calculates the response intensity as follows:

where the subscript "E1" means "Energy-type 1." The entire pulse energy can also be used as an intensity measure, such as in the mono-static and multi-static microwave imaging for breast cancer detection, and the intensity is given by:

where the subscript "E2" stands for "Energy-type 2." The peak value can be used as the response intensity measure due to the bipolar nature of the response at the focal point:

where the subscript "P" stands for "Peak," with P⁺ and P^' defined in Equations (18) and (19), respectively. Herein, the sign of the maximum amplitude is kept since the sign of the peak may also contain some information about the focal point. [0053] Peak-searching maximizes the output signal-to-noise ratio. An intuitive explanation is that, given the fact that the acoustic pulse is bipolar, if one assumes that the residual term e(t) is stationary, or its power is uniform over time, then the signal-to- noise ratio is maximized at the (positive or negative) peak of the acoustic pulse. In contrast, a conventional DAS fixes the samples to be summed up at the calculated arrival time. Due to phase distortions, the waveform at the calculated time may be far from the peak value.

[0054] The peak-to-peak difference can be employed as the response intensity for the focal point at location r :

/_PP = P⁺ -P- > O, (25)

where the subscript "PP" denotes the "Peak-to-Peak difference." Peak-to-peak difference has higher imaging contrast than peak value measure: the peak-to-peak difference of the bipolar pulse is approximately twice the absolute peak value; that is, the output signal power of the former is four times of the latter; yet the noise power of the former is only twice of the latter. The output Signal-to-Noise Ratio (SNR) can be doubled by using the peak-to-peak difference rather than the peak value. [0055] For example, the method of adaptive and robust detection for use in thermoacoustic tomography can include estimating a bipolar acoustic pulse generated by a focal point from one or more thermoacoustic signals using robust Capon beamforming, searching for two peaks of the bipolar acoustic pulse within an interval around a calculated arrival time ( the interval being determined from an estimated variance of arrival time differences, and the variance being generated by perturbing an array steering vector for allowing an uncertainty in phase measurement) and calculating a response intensity at the focal point by calculating a peak-to-peak difference of the two peaks. Referring again to FIG. 3, the image reconstruction unit 206 can enhance a contrast of the image by using the peak-to-peak difference for providing a response intensity of the focal point, wherein the waveform is bipolar having a positive peak corresponding to a compression pulse and a negative peak corresponding to a refraction pulse.

[0056] Where applicable, the present embodiments of the invention can be realized in hardware, software or a combination of hardware and software. Any kind of computer system or other apparatus adapted for carrying out the methods described herein are suitable. A typical combination of hardware and software can be a mobile communications device with a computer program that, when being loaded and executed, can control the mobile communications device such that it carries out the methods described herein. Portions of the present method and system may also be embedded in a computer program product, which comprises all the features enabling the implementation of the methods described herein and which when loaded in a computer system, is able to carry out these methods.

[0057] While the preferred embodiments of the invention have been illustrated and described herein, it will be clear that there can be other embodiments of the invention. Other embodiments will be apparent to those skilled in the art having the benefit of the description provided herein.

Claims

CLAIMSWhat is claimed is:

1. An adaptive and robust method of thermoacoustic tomography, comprising: receiving multiple thermoacoustic signals emitted from a focal point in response to radiation absorption of tissue induced by heating the tissue, wherein the at least one thermoacoustic signal includes a pulse component and a noise component; estimating a waveform of said pulse component by applying robust Capon beamforming to the multiple thermoacoustic signals for mitigating amplitude and phase distortion of the waveform caused by the noise component; and reconstructing a portion of an image of thermoacoustic response intensity corresponding to an absorption property of the tissue, the reconstruction based upon evaluating a response intensity from the waveform.

2. The method of claim 1 , wherein said estimating comprises: measuring an amplitude distortion from the multiple thermoacoustic signal; and introducing an uncertainity set in an array steering vector, thereby accommodating an uncertainty range during a phase measurement.

3. The method of claim 2, further comprising: searching a peak in the uncertainty range to estimate a phase distortion and comparing a location of the peak to an expected location; and providing a phase distortion correction to compensate the waveform for the phase distortion by shifting the waveform based on the phase distortion, wherein the phase distortion is associated with an acoustic speed variation corresponding to the tissue.

4. The method of claim 3, further comprising: enhancing a contrast of the image by using a peak-to-peak difference at the peak for providing response intensity of the focal point, wherein the waveform is bipolar having a positive peak corresponding to a compression pulse and a negative peak corresponding to a rarefaction pulse.

5. The method of claim 3, further comprising: scanning the focal point at multiple locations to cover an entire cross-section of the tissue for forming the image, wherein the scanning occurs at multiple heights.

6. The method of claim 3, further comprising: calculating an arrival time of the pulse component based on the peak; aligning the one or more thermoacoustic signals based on the arrival time; and compensating for a loss in amplitude due to a propagation decay associated with the tissue.

7. The method of claim 3, wherein the searching is performed within an interval around the arrival, and the interval is determined by a difference between a true arrival time and an actual arrival time, and the interval is based on an estimated variance of arrival time differences, and the interval is at least the duration of an expected pulse duration.

8. The method of claim 6, further comprising accounting for homogeneous acoustic properties of the tissue in the arrival time.

9. The method of claim 1, wherein a source of the radiation is one of microwave, laser, or electromagnetic stimulation.

10. An adaptive and robust method of thermoacoustic tomography, comprising: estimating an acoustic pulse generated by a focal point from one or more thermoacoustic signals using robust Capon beamforming; searching for a peak in the acoustic pulse for mitigating a phase aberration; and calculating a response intensity from the peak for forming one portion of an image.

11. The method of claim 10, wherein estimating a waveform includes: estimating an array steering vector based on covariance fitting; and obtaining an optimal beamformer weight vector from the estimated array steering vector to estimate a waveform of said pulse, wherein the array steering vector maximizes a power of the pulse subject to a dual constraint optimization based on the optimal beamformer weight vector.

12. The method of claim 10, wherein a first constraint of the dual optimization is a function of a covariance of the waveform, a power of the pulse in the waveform, and an array steering vector comprising element values for weighting the pulse.

13. The method of claim 10, wherein a second constraint of the dual optimization is an uncertainity set that is either a spherical uncertainty or an elliptical uncertainty, thereby accommodating an uncertainty range for the array steering vector during a phase measurement.

14. The method of claim 13, wherein an error of the uncertainity set is a user-provided parameter for adjusting a resolution and a contrast of the image, such that: reducing the uncertainity set provides an increased ability of the RCB to suppress a pulse interference for increasing a resolution of said image; and increasing the uncertainity set provides an increased robustness of the RCB to tolerate distortions due to a sample size of thermoacoustic signals.

15. The method of claim 13, further comprising: decreasing the uncertainity set for images having low resolution; and increasing the uncertainity set for images distorted by interference.

16. The method of claim 11 , wherein estimating an array steering vector further comprises: employing a Lagrange multiplier within the dual constraint optimization; applying an Eigenvalue decomposition to the covariance for creating a monotonically decreasing function of the Lagrange multiplier; and employing Newton's method to solve for the Lagrange multiplier.

17. The method of claim 16, further comprising eliminating a scaling ambiguity of the waveform by forcing a norm of the array steering vector to equal a number of employed transducers.

18. The method of claim 17, wherein estimating a waveform of the pulse comprises applying a weight vector to a pre-processed signal, and the weight vector is determined using the estimated array steering vector in a weight vector expression of a Capon beamformer, and the pre-processed signal is a gain-scaled and time-shifted version of the acoustic pulse.

19. The method of claim 13, wherein the uncertainity set is allowed for the array steering vector, thereby accommodating an uncertainty in phase measurement.

20. A method of adaptive and robust detection for use in thermoacoustic tomography comprising: estimating a bipolar acoustic pulse generated by a focal point from multiple thermoacoustic signals using robust Capon beamforming; searching for two peaks of the bipolar acoustic pulse within an interval around a calculated arrival time, wherein the interval is determined from an estimated variance of arrival time differences; and calculating a response intensity at the focal point by calculating a peak-to-peak difference of the two peaks.

21. An adaptive and robust imaging system for use in thermoacoustic tomography, comprising: a receiver for receiving one or more thermoacoustic signals emitted from a focal point in response radiation absorption of tissue induced by heating the tissue, wherein a thermoacoustic signal includes a pulse component and a noise component; a processor for estimating a waveform of said pulse component by applying robust Capon beamforming to said one or more thermoacoustic signals for mitigating an amplitude and phase distortion of said waveform due to said noise component; and an image reconstruction unit for evaluating a response intensity from said waveform for reconstructing a portion of an image of thermoacoustic response intensity corresponding to an absorption property of said tissue.

22. The system of claim 21 , wherein the processor further comprises: an analysis unit for measuring an amplitude distortion from said thermoacoustic signals, wherein said amplitude distortion is associated with a refraction property of said tissue.

23. The system of claim 22, wherein the processor further comprises: a peak detector for searching for a peak in said uncertainty range to estimate a phase distortion and comparing a location of said peak to an expected location; and a phase compensator for providing a phase distortion correction to compensate said waveform for said phase distortion by shifting said waveform by said phase distortion, wherein said phase distortion is associated with an acoustic speed variation due to said tissue.

24. The system of claim 23, wherein said image reconstruction unit enhances a contrast of said image by using a peak-to-peak difference at said peak for providing response intensity of said focal point, wherein said waveform is bipolar having a positive peak corresponding to a compression pulse and a negative peak corresponding to a rarefaction pulse.

25. The system of claim 23, wherein said receiver scans said focal point at one or more locations to cover an entire cross-section of said tissue for forming said image, wherein said scanning occurs at multiple heights.

26. The system of claim 23, wherein said processor calculates an arrival time of said pulse component based on said peak; aligns said one or more thermoacoustic signals based on said arrival time; and compensates for a loss in amplitude due to a propagation decay associated with said tissue.

27. The system of claim 26, wherein said peak detector searches for a peak within an interval around said arrival, and said interval is determined by a difference between a true arrival time and said arrival time

28. The system of claim 27, wherein said interval is based on an estimated variance of arrival time differences, and said interval is at least the duration of an expected pulse duration.