A High-Frequency Photoacoustic Computational Tomography Method Based on Sub-Nyquist Sampling

Abstract

This invention relates to the field of biomedical photoacoustic imaging technology and discloses a high-frequency photoacoustic computational tomography method based on sub-Nyquist sampling. The method includes: receiving a photoacoustic signal of a target via an ultrasonic transducer; passing the photoacoustic signal through a bandpass filter to remove noise and uninteresting spectra outside the spectrum of interest; sampling the photoacoustic signal after passing through the bandpass filter, flipping and copying the spectrum of interest in the high-frequency range to the low-frequency range to obtain a flipped spectrum; flipping and folding the flipped spectrum back into the high-frequency range to obtain a recovered spectrum, which is identical to the spectrum of interest; and performing photoacoustic computational tomography using the recovered spectrum to obtain an image of the target. This invention overcomes the limitations of existing PACT systems by introducing a high-frequency PACT method based on sub-Nyquist sampling, achieving high-frequency PACT with a relatively low sampling rate and obtaining high spatial resolution.

Description

Technical Field

[0001] This invention relates to the field of biomedical photoacoustic imaging technology, specifically to a high-frequency photoacoustic computational tomography method based on sub-Nyquist sampling. Background Technology

[0002] Photoacoustic computed tomography (PACT) is a rapidly developing biomedical imaging modality that has emerged in the last two decades. Based on the photoacoustic effect, PACT can provide high spatiotemporal resolution imaging of deep biological tissues, thus finding wide application in a range of biomedical fields, such as neurology, oncology, rheumatology, surgical navigation, and other medical imaging areas. High-resolution PACT can reveal more detailed structures within biological tissues. The spatial resolution of a PACT system is influenced by several system factors and is closely related to the characteristics of the ultrasound transducer (including center frequency, bandwidth, aperture size, and viewing angle). Spatial resolution can be easily improved by increasing the center frequency of the ultrasound transducer, thereby improving the image quality of PACT. According to the Nyquist sampling theorem, high-frequency PACT requires a data acquisition (DAQ) system with a high sampling rate analog-to-digital converter (ADC). In fact, the highest sampling rate supported by existing PACT systems is typically between 40 MHz and 60 MHz, resulting in a maximum supported Nyquist sampling frequency of no more than 30 MHz.

[0003] In existing technologies, DAQ systems with sampling rates higher than 60 MHz are typically used to achieve high-frequency PACT with bandwidths exceeding 30 MHz. Pan et al. used an oscilloscope with an 80 GHz sampling rate to sample photoacoustic signals detected by a micro-ring transducer array with a -6 dB bandwidth of 175 MHz. Vionnet et al. customized a DAQ system with a 125 MHz sampling rate to sample photoacoustic signals recorded by a 24 MHz transducer. Li et al. used a commercial imaging system with a 128 MHz sampling rate to sample photoacoustic signals received by a linear transducer array (bandwidth: 33 MHz) with a center frequency of 40 MHz. However, such high-sampling-rate DAQ systems are generally quite expensive. Inspired by previously reported delayed-triggered DAQ methods for high-frequency ultrasound imaging, Fu et al. proposed a new method called interleaved sampling photoacoustic imaging. This method enables high-frequency PACT at relatively low sampling rates, for example, sampling a 30 MHz transducer using a 41.67 MHz sampling rate. The proposed method requires two acquisitions, precisely offset by half the sampling period, both sampled at the same sampling rate (e.g., 41.67 MHz). These two sampled signals are interleaved to form a virtual sampled signal, equivalent to a signal sampled at twice the sampling rate (e.g., 83.33 MHz). Using this method, they achieved high-frequency photoacoustic imaging with an axial resolution of 63 μm and a lateral resolution of 91 μm. Summary of the Invention

[0004] To address the aforementioned technical problems, this invention provides a high-frequency photoacoustic computational tomography method based on sub-Nyquist sampling. By introducing a high-frequency PACT method based on sub-Nyquist sampling, the limitations of existing PACT systems are overcome, achieving high-frequency PACT with a relatively low sampling rate and obtaining higher spatial resolution.

[0005] To solve the above-mentioned technical problems, the present invention adopts the following technical solution:

[0006] A high-frequency photoacoustic computational tomography method based on sub-Nyquist sampling includes:

[0007] The photoacoustic signal from the target is received via an ultrasonic transducer; the center frequency of the photoacoustic signal is f. c The bandwidth is B;

[0008] The photoacoustic signal is passed through a bandpass filter to remove noise and uninteresting spectra outside the spectrum of interest.

[0009] Sub-Nyquist sampling is performed on the photoacoustic signal after passing through the bandpass filter, and the spectrum of interest is shifted to 0 to... The frequency range is used to obtain the shifted spectrum; where m is the number of shifts; and the sampling rate is... The shifted spectrum is then restored to the frequency range of the spectrum of interest to obtain the recovered spectrum; at this point, the equivalent sampling rate is... Furthermore, the recovered spectrum is exactly the same as the spectrum of interest. At this point, the sampling rate of the actual data acquisition system is... However, the recovered signal (recovered spectrum) is equivalent to using The signal is sampled at a sampling rate;

[0010] The image of the target is obtained by photoacoustic computational tomography through spectrum recovery.

[0011] Furthermore, the passband of the bandpass filter is arrive

[0012] Furthermore, the number of translations m = 1, and the signal bandwidth... Sampling rate The frequency range of the spectrum of interest is to Thus enabling the use The sub-Nyquist sampling rate can achieve the following results: The Nyquist sampling effect of the sampling rate.

[0013] Furthermore, in arrive The frequency range of interest and 0 to The frequency range shifted spectrum about Symmetry, when restoring the translated spectrum, the translated spectrum is along... Flip and fold back to The frequency range.

[0014] Furthermore, on a data acquisition system with a maximum sampling rate of 62.5MHz, sub-Nyquist sampling can be performed on a photoacoustic signal with a center frequency of 30MHz and a bandwidth of 20MHz at a sampling rate of 41.67MHz, but Nyquist sampling cannot be performed on the same signal at the same sampling rate of 62.5MHz. This is because the highest frequency of a photoacoustic signal with a center frequency of 30MHz and a bandwidth of 20MHz is 40MHz. Normal Nyquist sampling requires a sampling rate of 40*2 = 80MHz, while sub-Nyquist sampling only requires a sampling rate of 30*4 / 3 = 40MHz. Therefore, on a data acquisition system with a maximum sampling rate of 62.5MHz, sub-Nyquist sampling can be performed on this signal, but Nyquist sampling cannot.

[0015] Compared with the prior art, the beneficial technical effects of the present invention are:

[0016] This invention proposes an imaging technique based on sub-Nyquist sampling, enabling high-frequency PACT at a relatively low sampling rate, thereby reducing the need for expensive high-sampling-rate DAQ systems. Compared to traditional imaging methods, the imaging method of this invention exhibits a significant advantage in spatial resolution. To suppress ringing artifacts caused by this technique, this invention incorporates another signal sampling low-frequency components into the signal acquired by this sampling method. Under relatively low sampling rates (i.e., 41.67 MHz), this invention achieves an axial resolution of 22 μm using this imaging technique. It is understood that this resolution is the highest axial resolution ever achieved for PACT when the sampling rate does not exceed 60 MHz. Therefore, this invention provides a practically valuable solution for high-frequency PACT at limited sampling rates. Attached Figure Description

[0017] The above and other objects, features and advantages of the present invention will become more apparent from the following description of embodiments of the invention with reference to the accompanying drawings, in which:

[0018] Figure 1 The basic principle of the high-frequency PACT method based on sub-Nyquist sampling proposed in this invention is illustrated schematically. (a) is a flowchart of the invention; (b) is the original analog signal spectrum; (c) is the filtered analog signal spectrum, which is also the spectrum of interest; (d) is the sampled digital signal spectrum, where the high-frequency spectrum of interest is copied and flipped to the low-frequency range; and (e) is the high-frequency recovered spectrum obtained after flipping and recovery.

[0019] Figure 2 The performance of the imaging method in this invention is illustrated using numerical simulation. (a) is a schematic diagram of the simulation device, showing a microsphere with a diameter of 20 μm placed 8 mm away from the linear transducer array; (b) is the gold standard of the simulation obtained using sampling method 6 in Table 1; (c) to (g) are the simulation results obtained using sampling methods 1 to 5 in Table 1, respectively; (h) is the result of adding (c) and (d); (i) and (j) are comparisons of axial and lateral resolutions for different sampling methods.

[0020] Figure 3 The image results of the bionic experiment using body hair in this invention are illustrated schematically. (a) to (e) are the imaging results obtained using sampling methods 1 to 5 in Table 1, respectively; (f) is the result of adding (a) and (b).

[0021] Figure 4The image schematically illustrates the imaging results of a live experiment using a mouse ear, as described in this invention. (a) is a photograph of a mouse ear; (b) to (f) are the imaging results of blood vessels (indicated by white arrows) obtained using sampling methods 1 to 5 in Table 1; (g) is the sum of (b) and (c). The inset within the yellow dashed box is a magnified image of the blood vessels. (h) and (i) represent the blood vessel intensity distribution along the z-axis and x-axis, respectively, characterizing the full width at half maximum (FWHM). Detailed Implementation

[0022] A preferred embodiment of the present invention will now be described in detail with reference to the accompanying drawings.

[0023] Sub-Nyquist sampling, also known as bandpass sampling, allows bandpass signals to be sampled at a sampling rate less than twice that of the highest frequency component and is now widely used in ultrasound imaging to reduce computational complexity. The basic idea of ​​this invention is based on the fact that the frequency response of high-frequency ultrasound transducers is typically band-limited.

[0024] Suppose there is a center frequency f c A bandpass signal with bandwidth B can be sampled at f without aliasing, according to the sub-Nyquist sampling theorem. s The sampling rate is used to sample the bandpass signal. s The following conditions must be met:

[0025]

[0026] Where m is a definite expression satisfying f s Any positive integer greater than or equal to 2B.

[0027] Equation (1) can be rewritten as:

[0028]

[0029] The smaller the bandwidth, the larger the possible value of m, and the smaller the sampling rate f can be. s Since bandwidth is crucial for high-quality PACT, a smaller m should be chosen, such as m = 1 in this invention. In this case, the highest bandwidth supported by sub-Nyquist sampling is [missing value]. The corresponding f s yes The basic principle of the high-frequency PACT method based on sub-Nyquist sampling in this invention is as follows: Figure 1 As shown in (a), the sub-Nyquist sampling parameters used are: m = 1, First, a bandpass filter is used to filter out noise outside the range of interest and the spectrum of no interest, such as... Figure 1 (b) and Figure 1 As shown in (c) above. Next, the filtered photoacoustic signal (center frequency f) is...c bandwidth is Sampling is performed at a rate consistent with the highest frequency of the filtered photoacoustic signal.

[0030] After sampling, it is in the high-frequency range ( arrive The spectrum of interest will be inverted and copied to the low-frequency range (0 to 1) due to aliasing. ),like Figure 1 As shown in (d) in the diagram. Then, the spectrum to be copied to the low-frequency range is along... Flip and fold back to the high-frequency range, such as Figure 1 As shown in (e) in the diagram. At this point, the sampling rate has been doubled, i.e. Figure 1 In (e), the recovered spectrum and Figure 1 The spectrum of interest shown in (c) is exactly the same. Finally, the recovered signal is input into the algorithm for PACT image reconstruction.

[0031] Figure 1 This demonstrates that a signal with a bandwidth of 67% can be sampled at a sampling rate equal to the highest frequency of the signal and then accurately recovered through post-processing. The recovered signal is identical to the signal obtained by directly sampling the original signal at a sampling rate equal to twice the highest frequency of the signal (Nyquist sampling standard). The key to the sampling method in this invention is ensuring that the original analog signal (photoacoustic signal) is filtered by an analog bandpass filter before ADC sampling to eliminate noise and uninteresting spectra outside the spectrum of interest. Otherwise, this noise and uninteresting spectra will be aliased into the replicated spectrum after ADC sampling. In this case, it will be impossible to accurately recover the original spectrum of interest.

[0032] Table 1 shows the different sampling methods used for performance comparison in this invention.

[0033] Table 1

[0034]

[0035] Example 1: Numerical Simulation

[0036] Numerical simulations were performed by imaging microspheres with a diameter of 20 μm against a uniform background (velocity of sound: 1500 m / s). An analytical model was used to generate the photoacoustic signal from the microspheres. The generated photoacoustic signal was received by a linear transducer array with 256 elements and a 15 mm aperture. Each element had a center frequency of 30 MHz and a bandwidth of 20 MHz. The bandwidth characteristics of the transducer array were simulated by applying a Gaussian bandpass filter to the generated signal. The microspheres were placed 8 mm away from the transducer array. The received photoacoustic signal was first sampled at 1 GHz to be used as the raw signal in the numerical simulation. The signal was individually filtered using different filter configurations in Table 1 and downsampled to different sampling rates to compare the performance of this invention with other methods. In this numerical simulation, it was assumed that the highest sampling rate supported by DAQ was 62.5 MHz. Sampling rates of 83.33 MHz and 125 MHz were used to simulate interleaved sampling based on 41.67 MHz and 62.5 MHz, respectively. A first-order Butterworth high-pass filter and a third-order Butterworth low-pass filter are combined as a bandpass filter. A filtered back-projection algorithm is then used for PACT image reconstruction.

[0037] Numerical simulation devices such as Figure 2 As shown in (a), microspheres with a diameter of 20 μm are placed 8 mm away from the linear transducer array. Figure 2 (b) shows the imaging results of the microspheres obtained using sampling method 6 in Table 1. Although limited by transducer bandwidth and viewing angle, Figure 2 The reconstructed microspheres shown in (b) are not spherical, but Figure 2 (b) is still considered the gold standard in this numerical simulation. Because the non-circularity of the reconstructed microspheres is caused by hardware, the impact is consistent across all sampling methods and does not affect the comparison of the effects of different sampling methods here. Figure 2 Images (c) to (g) are the imaging results of the microspheres obtained using sampling methods 1 to 5 in Table 1. Figure 2 The negative artifacts shown in (b) to (g) are caused by the limited transducer bandwidth and viewing angle. Along Figure 2 The intensity distribution diagrams of the red dashed lines shown in (b) are respectively Figure 2 In the values ​​(i) and (j), the full width at half maximum (FWHM) values ​​represent the axial (z-axis) and lateral (x-axis) resolutions. It can be seen that... Figure 2 (c) is clearer than (e). This is because the signal obtained using sampling method 1 has higher frequency components than the signal obtained using sampling method 3. However, compared with... Figure 2 Compared to (b) in the middle, Figure 2(c) shows a more pronounced ringing artifact (indicated by the white arrow). This is because sampling mode 1 removes low-frequency signals, including only high-frequency signals. Figure 2 By fusing (d) and (c) in the above, a reconstruction result is obtained that simultaneously contains low-frequency (sampling mode 2) and high-frequency components (sampling mode 1), as shown below. Figure 2 As shown in (h) in the diagram. (and) Figure 2 Compared to (c) in the middle, Figure 2 (h) effectively suppresses ringing artifacts (indicated by the white arrow). It can be seen that the imaging method of this invention achieves the desired results (h). Figure 2 (c) in the middle and the results obtained based on interleaved sampling Figure 2 (f) and (g) are equally clear. Figure 2 The high similarity between (f) and (g) indicates that the spectrum of the microspheres is mainly within 35 MHz. These findings are strongly supported by both axial and lateral resolution. Figure 2 As shown in (i) and (j), the traditional sampling method (sampling method 3) yields the worst imaging results, with the widest microspheres. Numerical simulation results show that the high-frequency PACT method based on sub-Nyquist sampling proposed in this invention (sampling method 1) outperforms the traditional imaging method (sampling method 3) in spatial resolution and is comparable to the staggered sampling methods (sampling methods 4 and 5). To suppress ringing artifacts, it is necessary to add another sampling method (sampling method 2) that only acquires low-frequency components to the sampling method of this invention (sampling method 1).

[0038] Example 2: Phantom Experiment

[0039] To further verify the performance of the imaging method in this invention, firstly, a hair with a diameter of 10μm to 20μm was taken from the human forearm and embedded in an agar model. Then, the imaging method of this invention was used to image this model. The imaging system employed was a Vantage 256 system (Verasonics Inc., WA, USA), designed for real-time multi-channel (up to 256 channels) data acquisition and processing, supporting a maximum sampling rate of 62.5MHz. This system was connected to a commercial 256-channel linear transducer array (MS400, Verasonics Inc., Canada), with a center frequency of 30MHz, a bandwidth of 20MHz, and an aperture of approximately 15mm. The detected photoacoustic signals were individually filtered using various filter configurations and sampled at different rates (sampling methods 1 to 5 in Table 1) to compare the performance of the imaging method of this invention with other imaging methods. Other settings were the same as those in the numerical simulation. Figure 3 (a) to (e) in the table represent the hair imaging results obtained using different sampling methods (1 to 5 in Table 1), respectively. As expected, the imaging method of the present invention (see...) Figure 3 (a) provides a finer structure, but compared to traditional imaging methods (see [reference]). Figure 3 Compared to (c) in the previous example, there are more ringing artifacts (indicated by white arrows). Figure 3 (b) and (a) are combined. Figure 3 The ringing artifacts shown in (a) are well suppressed, as Figure 3 As shown in (f) in the figure. Compared to traditional imaging methods (see Figure 3 (c) in the middle, staggered sampling method (see Figure 3 (d) and (e) in the image also provide a clearer picture. Figure 3 The ringing artifact shown in (e) (indicated by the white arrow) is primarily caused by limited transducer bandwidth and viewing angle. It can be seen that the combined imaging method (see...) Figure 3 (f) performs better in terms of spatial resolution and image quality.

[0040] Example 3: In vivo experiment

[0041] To verify the performance of the imaging method of this invention, an in vivo experiment was conducted on the ears of nude mice. Figure 4 The mouse ear blood vessels shown in (a) were imaged using different sampling methods (sampling methods 1 to 5 in Table 1), and the results are as follows. Figure 4 As shown in (b) to (f) of the present invention. As expected, the imaging method of the present invention (see...) Figure 4 (b) in the middle provides higher spatial resolution, but compared with traditional imaging methods (see Figure 4 Compared to (d), ringing artifacts will appear. By... Figure 4 (c) and (b) are combined. Figure 4 The ringing artifacts shown in (b) are well suppressed, as Figure 4 As shown in (g) in the figure. The performance of the combined imaging method (see Figure 1). Figure 4 (g) in this context is comparable to the staggered sampling imaging method (see [link]). Figure 4 (e) and (f) in the table. These results were obtained. Figure 4 Further validation of the FWHM values ​​of the blood vessels shown in (h) and (i) in the figure.

[0042] In practice, the rising and falling edges of analog bandpass filters typically have finite slopes, meaning that noise and unwanted spectra outside the spectrum of interest cannot be completely eliminated. Therefore, the recovered spectrum will be contaminated to some extent. Figure 2 (h) and Figure 3The ringing artifact shown in (f) arises not only from the limited transducer bandwidth and viewing angle, but also from the finite slopes of the rising and falling edges of the analog bandpass filter. High-performance analog bandpass filters promise to reduce ringing artifacts.

[0043] It will be apparent to those skilled in the art that the present invention is not limited to the details of the exemplary embodiments described above, and that the invention can be implemented in other specific forms without departing from its spirit or essential characteristics. Therefore, the embodiments should be considered in all respects as exemplary and non-limiting, and the scope of the invention is defined by the appended claims rather than the foregoing description. Thus, all variations falling within the meaning and scope of equivalents of the claims are intended to be included within the present invention, and no reference numerals in the claims should be construed as limiting the scope of the claims.

[0044] Furthermore, it should be understood that although this specification describes embodiments, not every embodiment contains only one independent technical solution. This narrative style is merely for clarity. Those skilled in the art should consider the specification as a whole, and the technical solutions in each embodiment can also be appropriately combined to form other embodiments that can be understood by those skilled in the art.

Claims

1. A high-frequency photoacoustic computational tomography method based on sub-Nyquist sampling, comprising: The photoacoustic signal from the target is received via an ultrasonic transducer; the center frequency of the photoacoustic signal is f. c The bandwidth is B; The photoacoustic signal is passed through a bandpass filter to remove noise and uninteresting spectra outside the spectrum of interest. Sub-Nyquist sampling is performed on the photoacoustic signal after passing through the bandpass filter, and the spectrum of interest is shifted to 0 to... The frequency range is used to obtain the shifted spectrum; where m is the number of shifts; and the sampling rate is... The shifted spectrum is then restored back to the frequency range of the spectrum of interest to obtain the restored spectrum; at this point, the equivalent sampling rate is... Furthermore, the recovered spectrum is exactly the same as the spectrum of interest; Image of the target is obtained by photoacoustic computational tomography through spectrum recovery; The smaller the bandwidth, the larger the value of m that can be taken, and the larger the sampling rate f that can be taken. s The smaller the value; the more translations m=1, the greater the signal bandwidth. Sampling rate The frequency range of the spectrum of interest is to Thus enabling the use The sub-Nyquist sampling rate can achieve the following results: The Nyquist sampling effect of the sampling rate.

2. The high-frequency photoacoustic computational tomography method based on sub-Nyquist sampling according to claim 1, characterized in that, The passband of the bandpass filter is arrive 3. The high-frequency photoacoustic computational tomography method based on sub-Nyquist sampling according to claim 1, characterized in that, In arrive The frequency range of interest and 0 to The frequency range shifted spectrum about Symmetry, when restoring the translated spectrum, the translated spectrum is along... Flip and fold back to The frequency range.

4. The high-frequency photoacoustic computational tomography method based on sub-Nyquist sampling according to claim 1, characterized in that, On a data acquisition system with a maximum sampling rate of 62.5MHz, sub-Nyquist sampling can be performed on a photoacoustic signal with a center frequency of 30MHz and a bandwidth of 20MHz at a sampling rate of 41.67MHz, but Nyquist sampling cannot be performed on a photoacoustic signal with a center frequency of 30MHz and a bandwidth of 20MHz at a sampling rate of 62.5MHz.

Images