WO2016101382A1 - 微秒分辨空化时空分布的三维空化定量成像方法 - Google Patents

微秒分辨空化时空分布的三维空化定量成像方法 Download PDF

Info

Publication number
WO2016101382A1
WO2016101382A1 PCT/CN2015/071540 CN2015071540W WO2016101382A1 WO 2016101382 A1 WO2016101382 A1 WO 2016101382A1 CN 2015071540 W CN2015071540 W CN 2015071540W WO 2016101382 A1 WO2016101382 A1 WO 2016101382A1
Authority
WO
WIPO (PCT)
Prior art keywords
cavitation
dimensional
image
wide beam
quantitative
Prior art date
Application number
PCT/CN2015/071540
Other languages
English (en)
French (fr)
Inventor
万明习
丁婷
胡虹
杨淼
郭世放
Original Assignee
西安交通大学
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by 西安交通大学 filed Critical 西安交通大学
Priority to US15/115,682 priority Critical patent/US10332250B2/en
Publication of WO2016101382A1 publication Critical patent/WO2016101382A1/zh

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/0002—Inspection of images, e.g. flaw detection
    • G06T7/0012—Biomedical image inspection
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N29/00—Investigating or analysing materials by the use of ultrasonic, sonic or infrasonic waves; Visualisation of the interior of objects by transmitting ultrasonic or sonic waves through the object
    • G01N29/04—Analysing solids
    • G01N29/06—Visualisation of the interior, e.g. acoustic microscopy
    • G01N29/0654—Imaging
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N29/00—Investigating or analysing materials by the use of ultrasonic, sonic or infrasonic waves; Visualisation of the interior of objects by transmitting ultrasonic or sonic waves through the object
    • G01N29/02—Analysing fluids
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N29/00—Investigating or analysing materials by the use of ultrasonic, sonic or infrasonic waves; Visualisation of the interior of objects by transmitting ultrasonic or sonic waves through the object
    • G01N29/22—Details, e.g. general constructional or apparatus details
    • G01N29/26—Arrangements for orientation or scanning by relative movement of the head and the sensor
    • G01N29/265—Arrangements for orientation or scanning by relative movement of the head and the sensor by moving the sensor relative to a stationary material
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T11/00—2D [Two Dimensional] image generation
    • G06T11/003—Reconstruction from projections, e.g. tomography
    • G06T11/008—Specific post-processing after tomographic reconstruction, e.g. voxelisation, metal artifact correction
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00—Image enhancement or restoration
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/60—Analysis of geometric attributes
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N2291/00—Indexing codes associated with group G01N29/00
    • G01N2291/02—Indexing codes associated with the analysed material
    • G01N2291/024—Mixtures
    • G01N2291/02433—Gases in liquids, e.g. bubbles, foams
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N2291/00—Indexing codes associated with group G01N29/00
    • G01N2291/10—Number of transducers
    • G01N2291/106—Number of transducers one or more transducer arrays
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2200/00—Indexing scheme for image data processing or generation, in general
    • G06T2200/04—Indexing scheme for image data processing or generation, in general involving 3D image data
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/10—Image acquisition modality
    • G06T2207/10132—Ultrasound image
    • G06T2207/10136—3D ultrasound image
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/20—Special algorithmic details
    • G06T2207/20172—Image enhancement details
    • G06T2207/20182—Noise reduction or smoothing in the temporal domain; Spatio-temporal filtering

Definitions

  • the invention relates to the field of cavitation physics and application and ultrasonic imaging technology, which combines array Plane-by-plane wide beam cavitation detection, wide beam minimum variance adaptive beamforming and Nakagami parametric imaging algorithm to realize steady state free field and pulsation Microsecond time-spaced three-dimensional spatiotemporal distribution imaging and quantitative imaging of cavitation microbubble density under flow conditions.
  • Cavitation means that the cavitation nucleus in the liquid is activated by the applied energy (heat/force), and a series of kinetic processes such as oscillation, growth, shrinkage and even collapse of the tiny nucleus occur, which is the drug release in the biomedical field.
  • cavitation nucleation cavitation nucleation
  • cavitation bubble linear and nonlinear vibration cavitation bubble growth
  • cavitation rapid contraction to collapse rupture cavitation bubble dissipation
  • non-inertial cavitation Characteristic transient cavitation characterized by steady-state cavitation and inertial cavitation.
  • the cavitation initiation threshold the size of which depends on the static pressure of the medium liquid, the initial temperature, the structural state of the liquid itself, and the added cavitation nucleus in the liquid, so the liquid.
  • the generation of media hollowing has a certain randomness, but when the same environment and cavitation energy act, the shape and distribution of cavitation bubbles are reproducible.
  • effective cavitation detection and imaging methods need to be studied.
  • Optical detection imaging mainly includes high-speed/ultra-high-speed photography, sonoluminescence and sonochemiluminescence, which can observe the behavior of cavitation bubbles and spatio-temporal dynamic characteristics. It has the advantages of intuitiveness, good synchronization and high time resolution. On the one hand, the transparency of the medium is very high and it is not suitable for in-situ research. On the other hand, the obtained image is an overlap of information along the direction of light penetration.
  • the acoustic detection method is based on acoustic information generated during cavitation or cavitation microbubbles, including harmonics, subharmonics, superharmonics and broadband noise.
  • PCD Passive Cavitation Detection
  • ACD Active Cavitation Detection
  • PCI Passive Cavitation Imaging
  • Active Cavitation Imaging Active Cavitation Imaging
  • PCI obtains the two-dimensional spatial distribution of cavitation bubbles through the passive receiving of the array transducer and the reconstruction of the channel signal source.
  • the reconstruction algorithm is complex and the spatial resolution is not high.
  • ACI includes conventional B-ultrasound imaging and ultra-fast active cavitation imaging methods. Since the B-mode image is obtained by line-by-line scanning, there is a time difference between different scan lines of the same frame image, and the time resolution cannot reach the microsecond level.
  • Ultra-fast active cavitation imaging due to the emission of plane waves, has improved sensitivity and lateral resolution, and its time resolution is several hundred microseconds, which cannot meet the requirements of studying cavitation transient distribution.
  • cavitation bubbles need to be categorized, including cavitation quantification, cavitation size and density distribution.
  • the current cavitation quantification methods mainly include inertial cavitation dose and non-inertial cavitation dose, which is a relative measure of cavitation intensity by calculating the root mean square value of broadband noise or subharmonic amplitude in a specific frequency band.
  • the existing cavitation density detection method has a laser phase Doppler method, which mainly focuses on the distribution of cavitation bubbles in different sizes, and the density distribution of cavitation bubbles in different spatial positions has not been studied, and spatial information cannot be provided. .
  • the existing cavitation detection and imaging methods are limited to one-dimensional and two-dimensional, and the actual hollowedulized bubble distribution area is even larger throughout the focal region, and may exist in the acoustic wave propagation path in clinical applications such as focused ultrasound therapy. Other tissue media make the sound field distribution change and asymmetry occurs. Therefore, it is necessary to develop a quantitative imaging method for microsecond-resolved cavitation three-dimensional spatiotemporal distribution imaging and cavitation microbubble density.
  • there are relatively few studies on cavitation under flow conditions especially in pulsating flow conditions, and human blood flow is a pulsating flow, so it is necessary to study the three-dimensional spatiotemporal cavitation distribution under its conditions, especially in the pulsating flow cycle. Three-dimensional cavitation distribution at time.
  • Plane-by-plane wide beam detection cavitation is used to overcome the same frame cavitation image Disadvantages between the two, after each wide beam detection cavitation, the array transducer moves a unit position, and waits long enough for the medium cavitation nuclear distribution to return to the initial state before the cavitation energy source, and then Under the same cavitation energy excitation, the wide beam detection cavitation is again performed to gradually obtain a series of two-dimensional cavitation original RF data of different cell positions, and then combined with wide beam minimum variance adaptive beamforming (Synnevag, JF, A. Austeng) , et al. (2007).
  • Step 1 Using the array Plane-by-plane wide beam to detect cavitation, the cavitation is generated when the source energy temperature or pressure is continuously adjustable, and the cavitation signal is collected: the cavitation generating device includes the energy generating the energy field.
  • the detection device for the cavitation signal includes an array transducer that can transmit a wide beam and a parallel channel data acquisition and storage unit;
  • the synchronization signal generator generates a synchronization signal to respectively control the energy source device and
  • the array transducer the energy source device generates continuous variable energy excitation cavitation, the array transducer transmits a wide beam to detect cavitation, and the obtained cavitation echo signal is collected and stored by the parallel channel data acquisition and storage unit; Waiting for a long enough time to restore the media cavitation core distribution to the initial state, controlling the array transducer to move a cell position perpendicular to the array placement direction by a three-dimensional mechanical scanning device, and the cavitation energy source of the same parameter is again excited to generate cavitation
  • the synchronous array transducer transmits a wide beam to collect cavitation radio frequency data; repeating the above process, An energy source under cavitation conditions are the same parameters, a series of two-dimensional distribution image cavit
  • Step 2 Using wide beam minimum variance adaptive beamforming to process the original two-dimensional cavitation distribution image original RF data, and obtain a series of two-dimensional cavitation distribution images along different cell positions perpendicular to the array placement direction: first select Determining a certain target point in the two-dimensional cavitation imaging area of a certain unit position, calculating the effective aperture and delay according to the position of the target point to obtain the channel signal after receiving and focusing the target point; and then performing minimum variance adaptation on the channel signal Beam synthesis, realizing the optimal amplitude apodization of the target point channel signal and obtaining the optimal output; traversing all the target points of the two-dimensional cavitation imaging region to obtain the radio frequency data of the unit position two-dimensional cavitation imaging beam synthesis; the above process traversing A series of two-dimensional cavitation distribution radio frequency data is obtained for all cell locations. Finally, each two-dimensional cavitation radio frequency data is imaged by radio frequency imaging algorithm to obtain a series of high-resolution, high-signal-to-noise
  • Step 3 Perform a three-dimensional reconstruction on the one hand to obtain a three-dimensional cavitation distribution image on the one hand, and obtain a three-dimensional cavitation distribution image on the one hand, and obtain a beam obtained in step two on the other hand.
  • the synthesized radio frequency data is enveloped and then subjected to Nakagami parameter extraction to obtain a cavitation density distribution: using three-dimensional weight
  • the algorithm is used to perform three-dimensional display of two-dimensional cavitation distribution images of different unit positions.
  • the Nakagami parameter calculation is performed on the RF data after beam synthesis, and the two-dimensional cavitation density quantitative images of different unit positions are obtained, and then three-dimensional reconstruction is performed.
  • the algorithm obtains a three-dimensional display of cavitation density quantitative images.
  • Step 4 Change the parameters and the medium of the cavitation energy source device, and repeat steps 1 to 3 to obtain a cavitation and cavitation density quantitative three-dimensional image under different conditions: changing the cavitation energy source action time, and the time resolution can be reached. After a few microseconds, repeat steps one to three to obtain a cavitation and cavitation density quantitative three-dimensional sequence image that evolves over time; change the energy of the cavitation energy source, repeat steps one to three, and obtain cavitation and null with energy.
  • Densifying the quantitative three-dimensional sequence image changing the time delay between the excitation of the energy source device and the wide beam of the array transducer, repeating steps one to three, obtaining a cavitation and cavitation density quantitative three-dimensional sequence image dissipated over time; The medium for the action of the energy source is repeated steps one to three to obtain quantitative three-dimensional sequence images of cavitation and cavitation density under different media.
  • the Nakagami parameter extraction described in step 3 is as follows:
  • a three-dimensional cavitation imaging and cavitation density quantitative imaging method for spatiotemporal distribution of microsecond-resolved cavitation under pulsating flow conditions The principle is that the pulsation cycle time scale is much larger than the cavitation transient change time scale under pulsating flow conditions.
  • the cavitation kernel distribution can be restored to the initial state, and the Plane-by-plane wide beam detection cavitation is used to change the time delay of the pulsation pump and the cavitation energy source device, so that the pulsation period is not obtained.
  • the original RF data of a series of two-dimensional cavitation distribution images with different cell positions can be obtained by changing the position of the array transducer, and then combined with wide beam minimum variance adaptive beamforming, Nakagami parametric imaging and three-dimensional reconstruction algorithm to obtain pulsating flow.
  • a microsecond-resolved space-time distribution of a three-dimensional cavitation image and a quantitative image of cavitation microbubble density including the following steps:
  • Step 1 Simulate the blood vessel and the surrounding tissue by using the vascular phantom model.
  • the rubber tube is connected to both sides of the tube, and connected with the pulsation pump.
  • the vascular phantom and the three-dimensional mechanical scanning device are connected, and the cavitation is placed by three-dimensional movement.
  • the source acts on the area and causes cavitation to be created in the tube of the vascular phantom.
  • Step 2 Test the flow law under different pulse pump parameter settings: the pulsating liquid is a saline contrast agent microbubble solution, set the pulsation pump parameters, including the number of pulsation per minute and the flow rate, observe the Doppler spectrum and record the pulsation period and the law. Finally, the deionized water is flowed in the pipeline to rinse the physiological saline contrast agent microbubble solution in the pipeline.
  • Step 3 Using the array Plane-by-plane wide beam to detect cavitation, the cavitation is generated when the source energy temperature or pressure is continuously adjustable, and the cavitation signal is collected: setting the pulsation period and the flow rate, and editing the timing a pulsating flow device, a cavitation source device, and a cavitation signal detecting device, the pulsating flow device including a pulsation pump for generating a pulsating flow, a vascular phantom model, and a synchronous signal generator for controlling the operation thereof; the cavitation source device includes generating energy The energy source device of the field and the synchronous signal generator for controlling the timing, the detecting device for the cavitation signal comprises an array transducer for programmable transmission wide beam and a parallel channel data acquisition and storage unit; the synchronization signal generator generates a synchronization signal to respectively control the pulsation Pump, energy source device and array transducer, the pulsating pump causes the liquid to flow into the vascular phantom
  • Step 4 The position of the array transducer is unchanged, and the time delay between the trigger of the pulsation pump and the excitation of the cavitation energy source device is changed, and step 3 is repeated to obtain the original RF data of the sequence two-dimensional cavitation distribution image at each moment in the pulsation period:
  • the fixed array transducer changes the time delay of the pulsation pump and the energy source device according to the pulsation period and law recorded by the pulsating flow Doppler spectrum, and waits for a sufficiently long time to return the media cavitation core distribution to the initial state, repeating the steps Third, the original radio frequency data of the two-dimensional cavitation distribution image at different time points in the pulsation period is obtained.
  • Step 5 Change the position of the array transducer, repeat steps 3 and 4, and obtain the original RF data of the sequence two-dimensional cavitation distribution image at different spatial positions and their pulsation periods: control the array transducer along the three-dimensional mechanical scanning device Move a cell position perpendicular to the array placement direction, wait for a long enough time to return the media cavitation core distribution to the initial state, the cavitation energy source of the same parameter is again excited to generate cavitation, and the array transducer emits a wide beam acquisition cavitation
  • the radio frequency data repeating the above process, can obtain a series of two-dimensional cavitation distribution image original radio frequency data along different cell positions perpendicular to the array placement direction under the cavitation energy source condition of the same parameter.
  • the original radio frequency data of the two-dimensional cavitation distribution image at different time points in the pulsation period at the position can be obtained.
  • Step 6 The wide-beam minimum variance adaptive beamforming is used to process the acquired original RF data of the two-dimensional cavitation distribution image to obtain a series of two-dimensional cavitation distribution images along different cell positions perpendicular to the array placement direction: first select Determining a certain target point in the two-dimensional cavitation imaging area of a certain unit position, calculating the effective aperture and delay according to the position of the target point to obtain the channel signal after receiving and focusing the target point; and then performing minimum variance adaptation on the channel signal Beam synthesis, realizing the optimal amplitude apodization of the target point channel signal and obtaining the optimal output; traversing all the target points of the two-dimensional cavitation imaging region to obtain the radio frequency data of the unit position two-dimensional cavitation imaging beam synthesis; the above process traversing A series of two-dimensional cavitation distribution radio frequency data is obtained for all cell locations. Finally, each two-dimensional cavitation radio frequency data is imaged by radio frequency imaging algorithm to obtain a series of high-resolution, high-signal-to
  • Step 7 Perform a three-dimensional reconstruction on a series of two-dimensional cavitation distribution images obtained in step 6 along different unit positions perpendicular to the array placement direction to obtain a three-dimensional cavitation distribution image, and the beam obtained in step six on the other hand.
  • the synthesized radio frequency data is enveloped and then Nakagami parameter extraction is used to obtain the cavitation density distribution: the three-dimensional reconstruction algorithm is used to perform three-dimensional display of the two-dimensional cavitation distribution images of different unit positions; the radio frequency data after beam synthesis is first performed Nakagami
  • the parametric calculation results in a two-dimensional cavitation density quantitative image of different unit positions, and then a three-dimensional reconstruction algorithm is used to obtain a three-dimensional display of cavitation density quantitative images.
  • Step 8 Change the parameters and medium of the cavitation energy source device, repeat steps 3 to 7 to obtain cavitation and cavitation density quantitative three-dimensional images with different parameters under pulsating flow conditions: change the cavitation energy source action time, The time resolution can reach several microseconds, repeat steps three to seven, obtain the cavitation and cavitation density quantitative three-dimensional sequence image under the pulsating flow condition with the evolution of the cavitation source; change the energy of the cavitation energy source, repeat step three Up to seven, obtain a quantitative three-dimensional sequence image of cavitation and cavitation density under the condition of pulsating flow with energy; change the time delay between the excitation of the energy source device and the wide beam of the array transducer, repeat steps three to seven, and obtain Cavitation and cavitation density quantitative three-dimensional sequence images under time-pulsed pulsating flow conditions; changing the medium of cavitation energy source, repeat steps three to seven, and obtain quantitative three-dimensional sequence images of different media space-time and cavitation density under pulsating flow
  • the present invention first uses the array Plane-by-plane wide beam detection cavitation to acquire a series of two-dimensional spaces along different cell positions perpendicular to the array placement direction.
  • the original RF data of the image is distributed, and the two-dimensional cavitation distribution sequence image of different cell positions is obtained by combining the wide beam minimum variance adaptive beamforming, and then the Nakagami parameter algorithm is used to extract the quantitative image reflecting the cavitation density, and finally the three-dimensional reconstruction algorithm is used.
  • the time-series cavitation image can be obtained, and the change is produced.
  • the parameters and media of the cavitation energy source device can obtain cavitation and cavitation density quantitative three-dimensional sequence images, and the method has the potential to develop into a three-dimensional cavitation imaging standard method similar to sound field measurement.
  • the invention has the following advantages:
  • the present invention firstly uses Plane-by-plane wide beam detection cavitation, has microsecond time resolution, can observe cavitation bubble group transient distribution, and secondly combines wide beam minimum variance adaptive beamforming to obtain
  • the cavitation image has higher spatial resolution and signal-to-noise ratio; then the quantitative image of cavitation microbubble density distribution can be obtained by Nakagami parameter extraction; again, the two-dimensional cavitation distribution sequence of different cell positions is obtained by using 3D reconstruction algorithm.
  • the image is transformed into cavitation and cavitation density quantitative three-dimensional image; finally, changing the action time of the cavitation energy source device, the cavitation and cavitation density quantitative three-dimensional spatiotemporal distribution image with time evolution can be obtained, and the time resolution can reach several One microsecond; changing the time delay between excitation of the energy source device and the wide beam of the array transducer, obtaining cavitation and cavitation density quantitative three-dimensional spatiotemporal distribution images dissipated over time; changing the pulsation and excitation of the energy source device Time delay, the spatial and temporal distribution of cavitation and cavitation density at different time points in the pulsation period;
  • the energy of cavitation energy source can obtain cavitation and cavitation density quantitative three-dimensional sequence images with energy changes; change the media of cavitation energy source to obtain quantitative three-dimensional sequence images of different media space-time and cavitation density.
  • FIG. 1 is a schematic diagram of an array Plane-by-plane wide beam cavitation detecting apparatus of the present invention
  • Figure 3 is a flow chart of the Nakagami parameter imaging of the present invention.
  • Figure 4 is a schematic view of the pulsating flow microsecond cavitation imaging device of the present invention.
  • Figure 5 is a flow chart of the pulsating flow microsecond cavitation three-dimensional quantitative imaging of the present invention
  • FIG. 6, FIG. 7, and FIG. 8 are results of the present invention using a high-intensity focused ultrasound transducer as an energy source to excite tap water hollowing and cavitation density quantitative three-dimensional images.
  • the high-intensity focused ultrasound pulse length is 10 us-200 ms, and the electric power is 100 W.
  • Fig. 6(a) is a two-dimensional spatiotemporal distribution cavitation image
  • Fig. 6(b) is a two-dimensional space-time cavitation microbubble density quantitative image
  • Fig. 7 is an array transducer position from -5 mm to 6 mm when the pulse length is 20 ms. Spatial sequence cavitation imaging
  • Fig. 6(a) is a two-dimensional spatiotemporal distribution cavitation image
  • Fig. 6(b) is a two-dimensional space-time cavitation microbubble density quantitative image
  • Fig. 7 is an array transducer position from -5 mm to 6 mm when the pulse length is 20 ms
  • Fig. 8(a) is a three-dimensional cavitation image with a pulsation length of 20 ms
  • Fig. 8(b) is a three-dimensional cavitation microbubble density image with a pulsation length of 20 ms
  • FIG. 9 is a result of the time-spaced and cavitation density quantitative three-dimensional image of the pulsating flow excited by the high-intensity focused ultrasonic transducer as an energy source according to the present invention.
  • the high-intensity focused ultrasonic pulse has a length of 20 us and an electric power of 100 W, and the arrow indicates the direction of the ultrasonic wave.
  • cavitation energy source device 1 power amplifier 2
  • waveform generator 3 fully digitized Ultrasound device 4
  • first three-dimensional mechanical scanning device 5 second three-dimensional mechanical scanning device 5'
  • array transducer 6 sink 7
  • sound absorbing material 8 pulsating pump 9
  • storage pool 11 Waste liquid tank 12, latex tube 13.
  • an array Plane-by-plane wide beam cavitation detecting device includes a cavitation generating device and a cavitation signal detecting device.
  • the former includes a cavitation generating energy source device 1, a power amplifier 2, and a control timing synchronization.
  • a waveform generator 3 comprising a programmable transmit wide beam array transducer 6 of a fully digital ultrasound device 4 and a parallel channel data acquisition and storage unit. The first channel of the waveform generator 3 drives the power amplifier 2, and the cavitation energy source device generates an energy field.
  • the waveform of the waveform generator 3 can control the parameters of the energy source, including the energy action time and the working mode, and the time thereof.
  • the resolution can reach microseconds, and on the other hand, the size of the energy source can be controlled by the panel of the power amplifier 2.
  • the second channel of the waveform generator 3 drives the fully digital ultrasound device 4 to transmit a wide beam and acquire raw cavitation radio frequency data, the transmit and receive parameters of which can be implemented by programming the fully digital ultrasound device 4.
  • the waveform generator 3 realizes time synchronization between the cavitation generating means and the cavitation signal detecting means through the first channel and the second channel, and the time delay between the two can be set as required.
  • the energy source After the energy source generates cavitation and cavitation signal acquisition, waiting for a sufficiently long time to return the media cavitation core distribution to the initial state, moving a cell position in a direction perpendicular to the array transducer 6 by the three-dimensional mechanical scanning device 5.
  • the unit position can be set by the three-dimensional mechanical scanning device 5, the cavitation energy source device 1 of the same parameter is again excited to generate cavitation, and the synchronous array transducer 6 transmits a wide beam to collect cavitation radio frequency data, and the above process is repeated to obtain the same parameters.
  • a series of two-dimensional cavitation distribution image original RF data along different cell positions perpendicular to the array placement direction Under the condition of cavitation energy source, a series of two-dimensional cavitation distribution image original RF data along different cell positions perpendicular to the array placement direction. The above operation is carried out in the water tank 7, and the sound absorbing material 8 is placed on the bottom and side walls of the water tank 7.
  • the programmable fully digital ultrasound device 4 uses the array Plane-by-plane wide beam detection cavitation, transmits wide beam detection cavitation, obtains the original RF data of the cavitation signal, and then waits for a sufficient time to recover the media cavitation core distribution.
  • a unit position is moved by the three-dimensional mechanical scanning device 5 in a direction perpendicular to the array transducer 6, and the unit position size can be set by the three-dimensional mechanical scanning device 5, and the cavitation energy source device 1 of the same parameter is again excited.
  • the synchronous array transducer 6 transmits a wide beam to collect cavitation radio frequency data, repeats the above process, and obtains a series of two-dimensional cavitation distributions at different cell positions perpendicular to the array placement direction under the condition of cavitation energy source with the same parameters.
  • Image raw RF data
  • Nakagam parametric imaging algorithm is used to quantitatively image a series of two-dimensional cavitation radio frequency data, and the cavitation density parameters are extracted to obtain a series of two-dimensional cavitation density images along different cell positions perpendicular to the array placement direction;
  • Three-dimensional reconstruction algorithm is used to reconstruct the two-dimensional cavitation distribution image and the two-dimensional cavitation density image along different unit positions perpendicular to the array placement direction, and the three-dimensional cavitation distribution image and three-dimensional image under the action of cavitation energy source are obtained.
  • Cavitation density image
  • the probability density function of R obtained from the Nakagami statistical model is:
  • ⁇ ( ⁇ ) and U( ⁇ ) represent the gamma function and the unit step function, respectively, and m and ⁇ are the Nakagami parameters and scale parameters of the Nakagami statistical model, respectively;
  • a microsecond-resolved cavitation three-dimensional quantitative imaging apparatus for pulsating flow including a pulsating flow device, cavitation A generating device and a cavitation signal detecting device.
  • the pulsating flow device includes a pulsation pump 9, a vascular tissue phantom 10, and a waveform generator 3 that triggers its operation.
  • the cavitation generating means comprises a cavitation energy source device 1, a power amplifier 2 and a synchronous waveform generator 3 for controlling timing
  • the cavitation signal detecting device comprises a programmable transmit wide beam array transducer 6 of the fully digital ultrasonic device 4. And parallel channel data acquisition and storage units.
  • the first channel of the waveform generator 3 drives the power amplifier 2, and the cavitation energy source device generates an energy field.
  • the waveform of the waveform generator 3 can control the parameters of the energy source, including the energy action time and the working mode, and the time thereof.
  • the resolution can reach microseconds, and on the other hand, the size of the energy source can be controlled by the panel of the power amplifier 2.
  • the second channel of the waveform generator 3 drives the transmit wide beam of the fully digital ultrasound device 4 and acquires the original cavitation radio frequency data, the transmit and receive parameters of which can be implemented by programming the fully digitalized ultrasound device 4.
  • the Marker channel of the waveform generator 3 triggers the pulsation pump 9.
  • the waveform generator 3 realizes time synchronization between the pulsating flow device, the cavitation generating device, and the cavitation signal detecting device through the first channel, the second channel, and the Marker channel, and the time delay can be set according to requirements.
  • Both ends of the tube of the vascular tissue phantom 10 are connected to the latex tube 13 and connected to the pulsation pump 9, and the liquid is flowed from the reservoir 1 through the tube into the waste liquid pool 12, thereby simulating the pulsating blood flow.
  • the vascular tissue phantom 10 and the second three-dimensional mechanical scanning device 5' are connected, placed in the cavitation source action region by three-dimensional movement, and cavitation is generated in the tube of the vascular tissue phantom 10.
  • the all-digital ultrasonic device 4 uses the Doppler mode to test the Doppler spectrum under the parameters of the different pulsation pumps 9, and records the pulsation period and the law, thereby editing the synchronization waveform timing of the waveform generator 3.
  • the programmable transmit wide beam array transducer 6 of the all-digital ultrasound device 4 is coupled to the first three-dimensional mechanical scanning device 5, and after the energy source generates cavitation and cavitation signal acquisition, waits for a sufficiently long time for the media cavitation nucleus The distribution is restored to the initial state, and a unit position is moved by the first three-dimensional mechanical scanning device 5 in a direction perpendicular to the array transducer 6, and the unit position size can be set by the first three-dimensional mechanical scanning device 5, and the cavitation energy of the same parameter The source device 1 is again excited to generate cavitation, and the synchronous array transducer 6 transmits a wide beam to collect cavitation radio frequency data, and the above process is repeated to obtain a cavitation energy source with
  • the pulsating flow microsecond-resolved cavitation three-dimensional quantitative imaging process has the following steps:
  • vascular tissue phantom 10 Prepare a vascular tissue phantom 10, pour the gel solution into the mold, slowly remove the tube after cooling and solidification, and form a tube inside the solidified gel.
  • 3 Build a pulsating flow device, the tube of the vascular tissue phantom 10 is connected to the latex tube 13 on both sides, and is connected with the pulsation pump 9, and the liquid is flowed from the storage tank 11 into the waste liquid pool 12 through the pipeline to simulate the pulsating blood flow.
  • the liquid in the storage tank 11 is a saline contrast microbubble solution
  • the whole digital ultrasound Device 4 was changed to pulse Doppler mode to record the pulsation period and law under different pulsating pump parameters.
  • 5 Edit the channel 1, channel 2 and Marker waveforms of the waveform generator 3 according to the recorded pulsation period, respectively input to the power amplifier 2 to excite the cavitation source energy device 1 to generate cavitation, the all-digital ultrasound device 4 to control the array transducer 6
  • the transmitting and receiving plane wave acoustic signals and the pulsating pump 9 control the time delay of each time in the pulsating flow period and the cavitation energy source device 1.
  • the programmable fully digital ultrasound device 4 uses the array Plane-by-plane wide beam detection cavitation, transmits wide beam detection cavitation, obtains the original RF data of the cavitation signal, and then waits for a sufficiently long time to recover the media cavitation core distribution.
  • a unit position is moved by the three-dimensional mechanical scanning device 5 in a direction perpendicular to the array transducer 6, and the unit position size can be set by the three-dimensional mechanical scanning device 5, and the cavitation energy source device 1 of the same parameter is excited again.
  • the synchronous array transducer 6 transmits a wide beam to collect cavitation radio frequency data, repeats the above process, and obtains a series of two-dimensional cavitation distributions at different cell positions perpendicular to the array placement direction under the condition of cavitation energy source with the same parameters.
  • the minimum variance adaptive beamforming algorithm is used to beam synthesize the sequence two-dimensional cavitation original RF data, and obtain two-dimensional cavitation radio frequency data along different unit positions perpendicular to the array placement direction, followed by envelope detection , logarithmic compression and coordinate transformation to obtain a sequence two-dimensional cavitation distribution image; 8
  • the Nakagami parametric imaging algorithm is used to quantitatively image the two-dimensional cavitation radio frequency data, and the cavitation density parameter is extracted to obtain the two-dimensional cavitation density image of the sequence along different cell positions perpendicular to the array placement direction.
  • the three-dimensional reconstruction algorithm is used to The two-dimensional cavitation distribution image and the two-dimensional cavitation density image perpendicular to the cell placement direction are three-dimensionally reconstructed, and the three-dimensional cavitation distribution image and the three-dimensional cavitation density image under the cavitation energy source are obtained; The parameters and media of the cavitation energy source device are repeated in steps 5-9 to obtain a cavitation and cavitation density quantitative three-dimensional sequence image.
  • Example 1 Taking medium tap water as an example, the high-intensity focused ultrasound pulse length is 10 us-200 ms, and the electric power is 100 W.
  • the cavitation bubble group generated under the experimental conditions is cavitation and cavitation density quantitative three-dimensional imaging:
  • the tap water is injected into the water tank, the high-intensity focused ultrasonic transducer is fixed on the side wall of the water tank, and the wide beam is detected by the Sonix-Touch full digital ultrasonic device to detect cavitation microbubbles, and the array transducer is along the edge.
  • the high-intensity focused ultrasound is placed in the direction of propagation, and its center position corresponds to the focal region of the high-intensity focused ultrasound transducer, and the synchronization between the two is achieved by the two-channel arbitrary waveform generator AWG420; (2) through the waveform generator Set the high-intensity focused ultrasound (HIFU) pulse length to 10us-200ms, the power amplifier set the electric power to 100W, the waveform editing setting triggers the first channel of the power amplifier and the second channel that triggers the Sonix-Touch fully digital ultrasound device has a time delay of 1ms.
  • HIFU high-intensity focused ultrasound
  • the three-dimensional reconstruction algorithm is used to respectively sequence the two-dimensional cavitation distribution image with HIFU pulse length of 20ms.
  • the cavitation density distribution image is reconstructed in three dimensions to obtain a three-dimensional cavitation distribution image and a three-dimensional cavitation density image in the tap water.
  • FIG 7 shows the cross-sectional distribution of each cavitation bubble group with an array transducer position from -5 to 6 mm at a pulse length of 20 ms.
  • Figure 7 shows the spatial sequence cavitation distribution of the array transducer from -5 to 6 mm with a pulse length of 20 ms.
  • Figure 8 shows a three-dimensional cavitation image (a) and a three-dimensional cavitation microbubble density image (b) with a pulse length of 20 ms.
  • Example 2 Under pulsating flow conditions, the high-intensity focused ultrasound pulse length is 20 us, the electric power is 100 W, the pulsation frequency is 80 times/second, the flow rate is 5 cm/s, and the diameter of the pipe is 4 mm.
  • the cavitation bubble group generated under the experimental conditions is subjected to Cavitation and cavitation density quantitative three-dimensional imaging:
  • the pulsation frequency is set to 80 times/second
  • the flow rate is 5 cm/s
  • the storage solution is a saline contrast microbubble solution
  • the Sonix-Touch fully digital ultrasound device is set to pulse Doppler mode, and the pulsation parameter is recorded.
  • the Doppler spectrum is analyzed for its pulsation period and law to edit the waveform of the waveform generator AWG420.
  • the CH1, CH2 and Marker channels are input to the power amplifier, Sonix-Touch and pulsation pump respectively to achieve synchronization between the three.
  • the three-dimensional mechanical scanning device Multiscan 5800 is moved perpendicular to the array placement direction, and the moving unit distance is set to 1 mm, and the two-dimensional cavitation original radio frequency number is obtained along the sequence of obtaining different unit positions.
  • the envelope two-dimensional cavitation distribution image at each moment and different unit positions in the pulsation period is obtained by envelope detection, logarithmic compression and coordinate transformation.
  • Nakagami parametric imaging is adopted. The algorithm obtains a sequence two-dimensional cavitation density distribution image at each moment and different unit positions in the pulsation period;
  • Three-dimensional reconstruction algorithm is used to reconstruct the two-dimensional cavitation distribution image and the cavitation density distribution image respectively, and the three-dimensional cavitation distribution image and three-dimensional in the pipeline under high-intensity focused ultrasound pulse length 20us and electric power 100W are obtained. Cavitation density image.
  • Figure 9 gives the three-dimensional cavitation image (a) and the three-dimensional cavitation microbubble density image (b) under the conditions of Example 2, and it can be seen that the cavitation microbubbles are filled with the vascular tissue phantom tube and extend along both ends of the tube.
  • the array transducer moves a unit position, and the cavitation core distribution is restored. Then, under the same cavitation energy excitation, the wide beam detection cavitation is again obtained, and the space series of different unit positions are obtained. Cavitation of the original RF data, the time series can be obtained by changing the cavitation source energy action time, the time delay between the excitation of the energy source device and the wide beam emitted by the array transducer, and the time delay between the excitation of the pulsation pump and the energy source device.
  • Two-dimensional cavitation of original RF data and then combined with wide beam minimum variance adaptive beamforming, Nakagami parametric imaging, three-dimensional reconstruction algorithm to obtain quantitative images of microsecond-resolved cavitation three-dimensional spatiotemporal distribution images and cavitation microbubble density
  • the method combined with array Plane-by-plane wide beam detection cavitation, wide beam minimum variance adaptive beamforming and Nakagami parametric algorithm imaging not only have high spatial resolution, high signal to noise ratio imaging characteristics, but also can achieve steady state free field and pulsating flow conditions.

Landscapes

  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • General Health & Medical Sciences (AREA)
  • Health & Medical Sciences (AREA)
  • Pathology (AREA)
  • Biochemistry (AREA)
  • Analytical Chemistry (AREA)
  • Chemical & Material Sciences (AREA)
  • Immunology (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Theoretical Computer Science (AREA)
  • Acoustics & Sound (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Geometry (AREA)
  • Medical Informatics (AREA)
  • Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
  • Radiology & Medical Imaging (AREA)
  • Quality & Reliability (AREA)
  • Ultra Sonic Daignosis Equipment (AREA)
  • Investigating Or Analyzing Materials By The Use Of Ultrasonic Waves (AREA)

Abstract

一种微秒分辨空化时空分布的三维空化定量成像方法，使用宽波束检测空化后，阵列换能器移动一个单元位置，待空化核分布恢复，再在同样的空化能量激励下，再次宽波束检测空化，得到不同单元位置的空间系列二维空化原始射频数据，通过改变空化源能量作用时间、能量源装置激励与阵列换能器发射宽波束之间的时间延迟以及脉动泵与能量源装置激励之间的时间延迟，可得到时间序列二维空化原始射频数据，然后结合宽波束最小方差自适应波束合成、Nakagami参量成像、三维重建算法得到微秒分辨空化三维时空分布图像和空化微泡密度的定量图像，其时间分辨率可达到微秒。

Description

微秒分辨空化时空分布的三维空化定量成像方法 技术领域
本发明涉及空化物理与应用及超声成像技术领域，该方法结合阵列Plane-by-plane宽波束空化检测，宽波束最小方差自适应波束合成以及Nakagami参量成像算法，实现稳态自由场和脉动流条件下微秒分辨空化三维时空分布成像和空化微泡密度的定量成像。
背景技术
空化是指液体中的空化核在外加能量(热/力)的作用下被激活，出现微小泡核的振荡、生长、收缩乃至崩溃等一系列动力学过程，是生物医学领域中药物释放、基因转染、体外碎石、溶栓、止血、热疗以及肿瘤热消融等方面的主要机制。空化的过程包括以下几个阶段：空化成核、空化泡线性和非线性振动，空化泡生长、空化急速收缩至坍塌破裂以及空化泡消散，可分为以非惯性空化为特点的稳态空化和惯性空化为特点的瞬态空化。在液体介质中，产生空化的最小能量值称为空化起始阈值，其大小取决于媒介液体静态压、初始温度、液体本身的结构状态以及液体中外加的多样性空化核，因此液体媒介中空化的产生具有一定的随机性，但相同环境以及空化能量作用时，其空化泡群形状及分布具有可重复性。目前，为更好地研究不同媒介空化的产生机制以便更好地控制和利用空化，需要研究有效的空化检测与成像方法。
现有的空化检测与成像主要有光学和声学方法。光学检测成像主要包括通过高速/超高速摄影、声致发光以及声致化学发光等，可拍摄观察空化泡的行为以及时空动态特性，具有直观、同步性好、时间分辨率高的优点，缺点是一方面对媒介透光性要求很高且不适用于原位研究，另一方面所得图像是沿光穿透方向信息的重叠。声学检测方法是基于空化过程中或空化微泡产生的声信息，包括谐波、次谐波、超谐波和宽带噪声等，其中得到最广泛应用的就是被动空化检测(Passive Cavitation Detection,PCD)和主动空化检测(Active Cavitation Detection,ACD)。PCD利用换能器被动接收由空化微泡所产生的声散射信号，而ACD采用低压脉冲回波探测可能发生空化的区域，但PCD和ACD由于一般采用单阵元换能器，受限于有限的空间检测区域，无法提供空化微泡的空间分布。
在PCD和ACD基础上，使用二维阵列换能器作为空化检测换能器，发展出被动空化成像(Passive Cavitation Imaging,PCI)和主动空化成像(Active Cavitation Imaging, ACI)。由于超声空化具有瞬态特性，空化微泡的振动、坍塌破裂以及消散的时间都是微秒级，因此空化成像方法的时间分辨率需要达到微秒。同时，针对空化的瞬态特性，有必要得到空化微泡的时空分布，包括不同空化能量源作用时间和空化消散随时间的序列时空空化分布。PCI通过阵列换能器被动接收和通道信号源重建得到空化泡的二维空间分布，重建算法复杂且空间分辨率不高。ACI包括常规的B超成像和超快速主动空化成像方法。由于B超图像是通过逐线扫描得到的，同一帧图像不同扫描线之间存在时间差，且时间分辨率无法达到微秒级。而超快速主动空化成像由于发射的是平面波，其灵敏度以及横向分辨率有待改进，且其时间分辨率为几百个微秒，无法满足研究空化瞬态分布的要求。
在空化成像的基础上需要对空化泡进行定征，包括空化量化、空化尺寸及密度分布等。当前的空化量化方法主要有惯性空化剂量和非惯性空化剂量，是通过计算特定频段内宽带噪声或次谐波幅度的均方根值作为空化强度的一种相对度量，可分别衡量瞬态空化和稳态空化的相对大小，但一般针对PCD所采集到的一维射频数据，这种量化方法无法反应空化强度分布。现有的空化密度检测方法有激光相位多普勒法，该方法主要针对空化泡在不同尺寸上的分布，而关于空化泡在不同空间位置的密度分布还没有研究，无法提供空间信息。
目前已有的空化检测与成像方法局限于一维和二维，实际中空化泡的分布区域遍布整个焦域甚至更大，而且在临床应用中如聚焦超声治疗时，其声波传播路径上可能存在其他组织介质，使得声场分布发生变化而出现不对称性，因此有必要发展一种微秒分辨空化三维时空分布成像和空化微泡密度的定量成像方法。此外，对于流动条件尤其是脉动流条件下的空化研究相对较少，而人体的血流是一种脉动流，因此有必要研究其条件下的三维时空空化分布尤其是脉动流周期内不同时间点的三维空化分布。
发明内容
针对上述现有技术的缺陷以及微秒分辨空化时空分布的三维定量成像的必要性，本发明的目的在于提供一种稳态自由场和脉动流条件下具有微秒分辨空化时空分布的三维空化定量成像方法。
为了实现上述目的，本发明采取了如下的技术方案：
一种稳态自由场下微秒分辨空化时空分布的三维空化成像和空化密度定量成像方法，其原理是：
在稳态自由场条件下，采用Plane-by-plane宽波束检测空化，克服同帧空化图像空 间不同步的缺点，每次宽波束检测空化后，阵列换能器移动一个单元位置，并等待足够长的时间使得介质空化核分布恢复到空化能量源作用之前的初始状态，再在同样的空化能量激励下，再次宽波束检测空化，以此逐步得到不同单元位置的系列二维空化原始射频数据，再结合宽波束最小方差自适应波束合成(Synnevag,J.F.,A.Austeng,et al.(2007)."Adaptive beamforming applied to medical ultra sound imaging."Ieee Transactions on Ultrasonics Ferroelectrics and Frequency Control 54(8):1606-1613)、Nakagami参量成像和三维重建算法，得到稳态自由场下微秒分辨空化时空分布的三维空化图像和和空化微泡密度的定量图像。
步骤一、采用阵列Plane-by-plane宽波束检测空化，在源能量温度或压力连续可调的情况下激励空化的产生，采集空化信号：空化的产生装置包括产生能量场的能量源装置和控制时序的同步信号发生器；空化信号的检测装置包括可编程发射宽波束的阵列换能器和并行通道数据采集及存储单元；同步信号发生器产生同步信号分别控制能量源装置和阵列换能器，能量源装置产生连续可变能量激励空化的产生，阵列换能器发射宽波束对空化进行检测，得到的空化回波信号由并行通道数据采集及存储单元采集存储；等待一个足够长的时间使得媒介空化核分布恢复到初始状态，通过三维机械扫描装置控制阵列换能器沿垂直于阵列放置方向移动一个单元位置，相同参数的空化能量源再次激励产生空化，同步阵列换能器发射宽波束采集空化射频数据；重复以上过程，可得到相同参数的空化能量源条件下，沿垂直于阵列放置方向不同单元位置的一系列二维空化分布图像原始射频数据。
步骤二、采用宽波束最小方差自适应波束合成对采集到的二维空化分布图像原始射频数据进行处理，得到沿垂直于阵列放置方向不同单元位置的一系列二维空化分布图像：首先选定某一单元位置二维空化成像区域中的某个目标点，根据目标点的位置计算有效孔径及延时得到对目标点进行接收聚焦后的通道信号；然后对通道信号进行最小方差自适应波束合成，实现目标点通道信号的最优幅度变迹并得到最优输出；遍历该二维空化成像区域所有目标点得到该单元位置二维空化成像波束合成后的射频数据；以上过程遍历所有单元位置得到一系列二维空化分布的射频数据；最后通过射频成像算法对各二维空化射频数据进行成像，得到一系列高分辨、高信噪比的二维空化分布图像。
步骤三、对步骤二中得到的沿垂直于阵列放置方向不同单元位置的一系列二维空化分布图像一方面进行三维重建，得到三维空化分布图像，另一方面在步骤二中得到的波束合成后的射频数据取包络后进行Nakagami参量提取得到空化密度分布：利用三维重 建算法，将不同单元位置的序列二维空化分布图像进行三维显示；对波束合成后的射频数据先进行Nakagami参量计算，得到不同单元位置的序列二维空化密度定量图像，再经过三维重建算法，得到空化密度定量图像三维显示。
步骤四、改变产生空化的能量源装置的参数以及媒介，重复步骤一到三，得到不同条件下的空化以及空化密度定量三维图像：改变空化能量源作用时间，其时间分辨可达到几个微秒，重复步骤一到三，得到随时间演化的空化以及空化密度定量三维序列图像；改变空化能量源能量大小，重复步骤一到三，得到随能量变化的空化以及空化密度定量三维序列图像；改变能量源装置激励与阵列换能器发射宽波束之间的时间延迟，重复步骤一到三，得到随时间消散的空化以及空化密度定量三维序列图像；改变空化能量源作用的媒介，重复步骤一到三，得到不同媒介下空化以及空化密度定量三维序列图像。
步骤三中所述Nakagami参量提取，具体方法为：
1)将射频数据rf进行相关去噪处理，具体为：
a)选取一定尺寸大小的背景信号区域，计算该区域平均能量P；
b)对rf分别叠加平均能量为P的随机高斯白噪声n1和n2，得到S1和S2；
c)计算S1和S2的相关系数，并给定阈值Th，对相关系数进行阈值化处理后与rf加权，得到rfdenoise如下，其中corrcoef(S1,S2)表示信号相关系数：
Figure PCTCN2015071540-appb-000001
2)将去噪后的空化射频数据rfdenoise进行Hilbert解调，得到包络信号，记为R；
3)计算Nakagami参量如下：
Figure PCTCN2015071540-appb-000002
其中E(·)代表统计均值。
一种脉动流条件下的微秒分辨空化时空分布的三维空化成像和空化密度定量成像方法，其原理是在脉动流条件下，由于脉动周期时间尺度远大于空化瞬态改变时间尺度，通过脉动周期倍数的等待时间可使得空化核分布恢复到初始状态，采用Plane-by-plane宽波束检测空化，改变脉动泵和空化能量源装置的时间延迟，可得到脉动周期内不同时 刻，而改变阵列换能器位置，可得到不同单元位置的系列二维空化分布图像原始射频数据，然后再结合宽波束最小方差自适应波束合成、Nakagami参量成像和三维重建算法，得到脉动流下微秒分辨空化时空分布的三维空化图像和空化微泡密度的定量图像，包括以下步骤：
步骤一、利用仿血管组织体模模拟血管及其周围的组织，管道两边连接橡皮管，并与脉动泵相连，连接仿血管组织体模与三维机械扫描装置，通过三维移动将其放置于空化源作用区域，并使空化产生在仿血管组织体模的管道内。
步骤二、测试不同脉动泵参数设置下的流动规律：脉动液体为生理盐水造影剂微泡溶液，设置脉动泵参数，包括每分钟脉动次数以及流量，观察多普勒谱并记录脉动周期和规律，最后管道内流动去离子水将管道内的生理盐水造影剂微泡溶液冲洗干净。
步骤三、采用阵列Plane-by-plane宽波束检测空化，在源能量温度或压力连续可调的情况下激励空化的产生，采集空化信号：设置脉动周期和流量，并以此编辑时序，控制脉动流装置、空化源装置和空化信号检测装置，脉动流装置包括产生脉动流的脉动泵、仿血管组织体模和控制其工作的同步信号发生器；空化源装置包括产生能量场的能量源装置和控制时序的同步信号发生器，空化信号的检测装置包括可编程发射宽波束的阵列换能器和并行通道数据采集及存储单元；同步信号发生器产生同步信号分别控制脉动泵、能量源装置和阵列换能器，脉动泵使液体按照一定脉动周期和流量流入流出仿血管组织体模，能量源装置产生连续可变能量激励空化的产生，阵列换能器发射宽波束对空化进行检测，得到的空化回波信号由并行通道数据采集及存储单元采集存储。
步骤四、阵列换能器位置不变，改变脉动泵触发与空化能量源装置激励之间的时间延迟，重复步骤三，得到脉动周期内各时刻的序列二维空化分布图像原始射频数据：固定阵列换能器，依据脉动流多普勒谱记录的脉动周期和规律，改变脉动泵与能量源装置的时间延迟，等待一个足够长的时间使得媒介空化核分布恢复到初始状态，重复步骤三，得到脉动周期内不同时间点的二维空化分布图像原始射频数据。
步骤五、改变阵列换能器位置，重复步骤三、四，得到不同空间位置及其脉动周期内各时刻的序列二维空化分布图像原始射频数据：通过三维机械扫描装置控制阵列换能器沿垂直于阵列放置方向移动一个单元位置，等待一个足够长的时间使得媒介空化核分布恢复到初始状态，相同参数的空化能量源再次激励产生空化，阵列换能器发射宽波束采集空化射频数据，重复以上过程，可得到相同参数的空化能量源条件下，沿垂直于阵列放置方向不同单元位置的一系列二维空化分布图像原始射频数据。不同单元位置时， 重复步骤四，可得到该位置下脉动周期内不同时间点的二维空化分布图像原始射频数据。
步骤六、采用宽波束最小方差自适应波束合成对采集到的二维空化分布图像原始射频数据进行处理，得到沿垂直于阵列放置方向不同单元位置的一系列二维空化分布图像：首先选定某一单元位置二维空化成像区域中的某个目标点，根据目标点的位置计算有效孔径及延时得到对目标点进行接收聚焦后的通道信号；然后对通道信号进行最小方差自适应波束合成，实现目标点通道信号的最优幅度变迹并得到最优输出；遍历该二维空化成像区域所有目标点得到该单元位置二维空化成像波束合成后的射频数据；以上过程遍历所有单元位置得到一系列二维空化分布的射频数据；最后通过射频成像算法对各二维空化射频数据进行成像，得到一系列高分辨、高信噪比的二维空化分布图像。
步骤七、对步骤六中得到的沿垂直于阵列放置方向不同单元位置的一系列二维空化分布图像一方面进行三维重建，得到三维空化分布图像，另一方面在步骤六中得到的波束合成后的射频数据取包络后进行Nakagami参量提取得到空化密度分布：利用三维重建算法，将不同单元位置的序列二维空化分布图像进行三维显示；对波束合成后的射频数据先进行Nakagami参量计算，得到不同单元位置的序列二维空化密度定量图像，再经过三维重建算法，得到空化密度定量图像三维显示。
步骤八、改变产生空化的能量源装置的参数以及媒介，重复步骤三到七，得到脉动流条件下不同参数时的空化以及空化密度定量三维图像：改变空化能量源作用时间，其时间分辨可达到几个微秒，重复步骤三到七，得到随空化源作用时间演化的脉动流条件下空化以及空化密度定量三维序列图像；改变空化能量源能量大小，重复步骤三到七，得到随能量变化的脉动流条件下空化以及空化密度定量三维序列图像；改变能量源装置激励与阵列换能器发射宽波束之间的时间延迟，重复步骤三到七，得到随时间消散的脉动流条件下空化以及空化密度定量三维序列图像；改变空化能量源作用的媒介，重复步骤三到七，得到脉动流条件下不同媒介时空化以及空化密度定量三维序列图像。
针对稳态自由场和脉动流条件下的空化时空分布研究，本发明首先使用阵列Plane-by-plane宽波束检测空化，采集得到沿垂直于阵列放置方向不同单元位置的一系列二维空化分布图像原始射频数据，再结合宽波束最小方差自适应波束合成得到不同单元位置的二维空化分布序列图像，然后使用Nakagami参量算法提取得到反映空化密度的定量图像，最后使用三维重建算法得到微秒分辨空化的三维定量图像，改变空化源能量作用时间以及其与宽波束检测之间的时间延迟，可得到时间序列空化图像，而改变产 生空化的能量源装置的参数以及媒介，可得到空化以及空化密度定量三维序列图像，该方法有潜力发展成为类似于声场测量的三维空化成像标准方法。
本发明与现有技术相比，具有下列优点：
与现有技术相比，本发明首先采用Plane-by-plane宽波束检测空化，具有微秒时间分辨率，可观察空化泡群瞬态分布；其次结合宽波束最小方差自适应波束合成得到的空化图像具有较高的空间分辨率和信噪比；然后通过Nakagami参量提取可得到空化微泡密度分布的定量图像；再次，使用三维重建算法将不同单元位置的二维空化分布序列图像转化为空化以及空化密度定量三维图像；最后，改变产生空化的能量源装置作用时间，可得到随时间演化的空化以及空化密度定量三维时空分布图像，其时间分辨可达到几个微秒；改变能量源装置激励与阵列换能器发射宽波束之间的时间延迟，可得到随时间消散的空化以及空化密度定量三维时空分布图像；改变脉动与能量源装置激励之间的时间延迟，可得到脉动周期内不同时间点的空化及空化密度定量时空分布图像；改变空化能量源能量大小，可得到随能量变化的空化以及空化密度定量三维序列图像；改变空化能量源作用的媒介，得到不同媒介时空化以及空化密度定量三维序列图像。
附图说明
图1是本发明阵列Plane-by-plane宽波束空化检测装置示意图；
图2是本发明的微秒分辨空化的三维定量成像流程图；
图3是本发明的Nakagami参量成像流程图；
图4是本发明脉动流微秒空化成像装置示意图；
图5是本发明的脉动流微秒空化三维定量成像流程图；
图6、图7、图8是本发明以高强度聚焦超声换能器作为能量源激励自来水中空化以及空化密度定量三维图像结果，高强度聚焦超声脉冲长度10us-200ms，电功率100W，箭头所指为超声波作用方向。图6(a)为二维时空分布空化图像，图6(b)为二维时空空化微泡密度定量图像，图7是脉冲长度为20ms时阵列换能器位置从-5mm到6mm的空间序列空化成像，图8(a)是脉动长度为20ms时三维空化图像，图8(b)是脉动长度为20ms时三维空化微泡密度图像；
图9是本发明以高强度聚焦超声换能器作为能量源激励脉动流时空化以及空化密度定量三维图像结果，高强度聚焦超声脉冲长度20us，电功率100W，箭头所指为超声波作用方向。
图中所示的序号为：空化能量源装置1，功率放大器2，波形发生器3，全数字化 超声设备4，第一三维机械扫描装置5，第二三维机械扫描装置5'，阵列换能器6，水槽7，吸声材料8，脉动泵9，仿血管组织体模10，储存池11，废液池12，乳胶管13。
具体实施方式
以下结合附图和实施例对本发明作详细说明。
参见图1，一种阵列Plane-by-plane宽波束空化检测装置，包括空化产生装置和空化信号检测装置，前者包括产生空化的能量源装置1、功率放大器2和控制时序的同步波形发生器3，后者包括全数字化超声设备4的可编程发射宽波束阵列换能器6和并行通道数据采集及存储单元。波形发生器3的第一通道驱动功率放大器2，激励空化能量源装置产生能量场，一方面通过波形发生器3的波形编辑可控制能量源的参数，包括能量作用时间以及工作模式，其时间分辨率可达到微秒，另一方面通过功率放大器2的面板可控制能量源的大小。波形发生器3的第二通道驱动全数字化超声设备4发射宽波束和采集原始空化射频数据，其发射接收参数可通过对全数字化超声设备4编程实现。波形发生器3通过第一通道和第二通道实现空化产生装置和空化信号检测装置之间的时间同步，可根据需求设置二者之间的时间延迟。当能量源产生空化和空化信号采集之后，等待一个足够长的时间使得媒介空化核分布恢复到初始状态，通过三维机械扫描装置5沿垂直于阵列换能器6的方向移动一个单元位置，单元位置可通过三维机械扫描装置5进行设置，相同参数的空化能量源装置1再次激励产生空化，同步阵列换能器6发射宽波束采集空化射频数据，重复以上过程，得到相同参数的空化能量源条件下，沿垂直于阵列放置方向不同单元位置的一系列二维空化分布图像原始射频数据。以上操作在水槽7中进行，并在水槽7的底部和侧壁放置吸声材料8。
参见图2，微秒分辨空化的三维定量成像流程，其步骤如下：
①采用阵列Plane-by-plane宽波束检测空化，可编程全数字化超声设备4发射宽波束检测空化，得到空化信号原始射频数据，接着等待一个足够长的时间使得媒介空化核分布恢复到初始状态，通过三维机械扫描装置5沿垂直于阵列换能器6的方向移动一个单元位置，单元位置大小可通过三维机械扫描装置5进行设置，相同参数的空化能量源装置1再次激励产生空化，同步阵列换能器6发射宽波束采集空化射频数据，重复以上过程，得到相同参数的空化能量源条件下，沿垂直于阵列放置方向不同单元位置的一系列二维空化分布图像原始射频数据；
②采用最小方差自适应波束合成算法对所述二维空化原始射频数据进行波束合成，得到沿垂直于阵列放置方向不同单元位置的一系列二维空化射频数据，随后包络检波、 对数压缩以及坐标变换得到一系列二维空化分布图像；
③采用Nakagam参量成像算法对一系列二维空化射频数据进行定量成像，提取空化密度参量，得到沿垂直于阵列放置方向不同单元位置的一系列二维空化密度图像；
④采用三维重建算法，对沿垂直于阵列放置方向不同单元位置的序列二维空化分布图像和二维空化密度图像进行三维重建，得到空化能量源作用下的三维空化分布图像和三维空化密度图像；
⑤改变产生空化的能量源装置的参数以及媒介，重复步骤①-④，得到空化以及空化密度定量三维序列图像。
参见图3，Nakagami参量成像算法流程，其步骤如下：
①取宽波束最小方差自适应波束合成后的射频数据记为rf，计算背景区域平均能量，记为P；
②对rf叠加平均能量均为P的高斯白噪声信号得到加噪信号S1和S2；
③给定阈值Th，对S1和S2进行相关性分析并对相关系数阈值化处理得到二值化相关系数矩阵Coefcorr；
④将rf与③中得到的二值化矩阵Coefcorr进行加权得到相关去噪后的射频数据rfdenoise：
Figure PCTCN2015071540-appb-000003
⑤对④中得到的rfdenoise进行Hilbert变化得到包络信号R；
由Nakagami统计模型得到R的概率密度函数为：
Figure PCTCN2015071540-appb-000004
其中Γ(·)和U(·)分别代表gamma函数和单位阶跃函数，m和Ω分别为Nakagami统计模型的Nakagami参量和尺度参量；
⑥利用公式
Figure PCTCN2015071540-appb-000005
计算Nakagami参量。
参见图4，一种脉动流的微秒分辨空化三维定量成像装置，包括脉动流装置、空化 产生装置和空化信号检测装置。脉动流装置包括脉动泵9、仿血管组织体模10和触发其工作的波形发生器3。空化产生装置包括产生空化的能量源装置1、功率放大器2和控制时序的同步波形发生器3，而空化信号检测装置包括全数字化超声设备4的可编程发射宽波束阵列换能器6和并行通道数据采集及存储单元。波形发生器3的第一通道驱动功率放大器2，激励空化能量源装置产生能量场，一方面通过波形发生器3的波形编辑可控制能量源的参数，包括能量作用时间以及工作模式，其时间分辨率可达到微秒，另一方面通过功率放大器2的面板可控制能量源的大小。波形发生器3的第二通道驱动全数字化超声设备4的发射宽波束和采集原始空化射频数据，其发射接收参数可通过对全数字化超声设备4编程实现。波形发生器3的Marker通道触发脉动泵9。波形发生器3通过第一通道、第二通道和Marker通道实现脉动流装置、空化产生装置和空化信号检测装置之间的时间同步，可根据需求设置时间延迟。将仿血管组织体模10的管道两端连接乳胶管13，并与脉动泵9相连，将液体从储存池11经管道流入废液池12，以此模拟脉动血流。连接仿血管组织体模10与第二三维机械扫描装置5'，通过三维移动将其放置于空化源作用区域，并使空化产生在仿血管组织体模10的管道内。全数字化超声设备4采用多普勒模式，测试不同脉动泵9参数下的多普勒谱，记录脉动周期和规律，以此编辑波形发生器3的同步波形时序。全数字化超声设备4的可编程发射宽波束阵列换能器6与第一三维机械扫描装置5相连，当能量源产生空化和空化信号采集之后，等待一个足够长的时间使得媒介空化核分布恢复到初始状态，通过第一三维机械扫描装置5沿垂直于阵列换能器6的方向移动一个单元位置，单元位置大小可通过第一三维机械扫描装置5进行设置，相同参数的空化能量源装置1再次激励产生空化，同步阵列换能器6发射宽波束采集空化射频数据，重复以上过程，得到相同参数的空化能量源条件下，沿垂直于阵列放置方向不同单元位置的一系列二维空化分布图像原始射频数据。以上操作在水槽7中进行，并在水槽7的底部和侧壁放置吸声材料8。
参见图5，脉动流微秒分辨空化三维定量成像流程，其步骤如下：
①制备仿血管组织体模10，将凝胶溶液倒入模具，待冷却凝固后慢慢取出管子，在凝固凝胶内部形成管道。②定位仿血管组织体模10，先连接仿血管组织体模10和第二三维机械扫描装置5'，通过三维移动将其放置于空化源作用区域，并使空化产生在仿血管组织体模10的管道内。③搭建脉动流装置，仿血管组织体模10的管道两边连接乳胶管13，并与脉动泵9相连，从储存池11中将液体经管道流入废液池12，模拟脉动血流。④测量脉动流流动，储存池11中液体为生理盐水造影微泡溶液，将全数字化超声 设备4改为脉冲多普勒方式，记录不同脉动泵参数下的脉动周期和规律。⑤根据记录的脉动周期，编辑波形发生器3的通道1、通道2和Marker波形，分别输入给功率放大器2激励空化源能量装置1产生空化、全数字化超声设备4控制阵列换能器6发射接收平面波声信号、以及脉动泵9控制脉动流周期内各时刻与空化能量源装置1的时间延迟。⑥采用阵列Plane-by-plane宽波束检测空化，可编程全数字化超声设备4发射宽波束检测空化，得到空化信号原始射频数据，接着等待一个足够长的时间使得媒介空化核分布恢复到初始状态，通过三维机械扫描装置5沿垂直于阵列换能器6的方向移动一个单元位置，单位位置大小可通过三维机械扫描装置5进行设置，相同参数的空化能量源装置1再次激励产生空化，同步阵列换能器6发射宽波束采集空化射频数据，重复以上过程，得到相同参数的空化能量源条件下，沿垂直于阵列放置方向不同单元位置的一系列二维空化分布图像原始射频数据；⑦采用最小方差自适应波束合成算法对序列二维空化原始射频数据进行波束合成，得到沿垂直于阵列放置方向不同单元位置的序列二维空化射频数据，随后包络检波、对数压缩以及坐标变换得到序列二维空化分布图像；⑧采用Nakagami参量成像算法对序列二维空化射频数据进行定量成像，提取空化密度参量，得到沿垂直于阵列放置方向不同单元位置的序列二维空化密度图像；⑨采用三维重建算法，对沿垂直于阵列放置方向不同单元位置的序列二维空化分布图像和二维空化密度图像进行三维重建，得到空化能量源作用下的三维空化分布图像和三维空化密度图像；⑩改变产生空化的能量源装置的参数以及媒介，重复步骤⑤-⑨，得到空化以及空化密度定量三维序列图像。
实施例1.以介质自来水为例,高强度聚焦超声脉冲长度10us-200ms，电功率100W，对该实验条件下产生的空化泡群进行空化以及空化密度定量三维成像：
(1)采用图1的实验装置，水槽中注入自来水，高强度聚焦超声换能器固定于水槽侧壁，采用Sonix-Touch全数字化超声设备发射宽波束探测空化微泡，阵列换能器沿着高强度聚焦超声传播方向放置，且使得其中心位置对应于高强度聚焦超声换能器的焦域，通过双通道任意波形发生器AWG420实现两者之间的同步；(2)通过波形发生器设置高强度聚焦超声(HIFU)脉冲长度为10us-200ms，功率放大器设置电功率为100W，波形编辑设置触发功率放大器的第一通道和触发Sonix-Touch全数字化超声设备的第二通道的时间延迟为1ms；(3)采用阵列Plane-by-plane宽波束检测空化，通过三维机械扫描装置Multiscan 5800沿垂直于阵列放置方向移动，移动单元设置为0.5mm，得到不同单元位置的序列二维空化原始射频数据；(4)采用最小方差自适应波束合成算法对序 列二维空化原始射频数据进行波束合成，得到不同单元位置的序列二维空化射频数据；(5)针对序列二维空化射频数据，一方面采用包络检波、对数压缩以及坐标变换得到序列二维空化分布图像，一方面采用Nakagami参量成像算法得到序列二维空化密度分布图像；(6)采用三维重建算法，分别对HIFU脉冲长度为20ms的序列二维空化分布图像和空化密度分布图像进行三维重建，得到自来水中三维空化分布图像和三维空化密度图像。
当阵列换能器位置为0mm时，HIFU脉冲长度为10us-200ms的空化时空和空化密度时空分布结果参见图6，根据分析，空化微泡出现在焦域，这是因为焦域声压最高，但泡群中心位置并不是严格位于0mm焦点位置，而是偏向焦后区域，这是因为声场产生的辐射力会对空化微泡产生作用力，使其沿声传播的方向发生移动，当脉冲长度为5ms时空化泡群逐渐在焦后形成枝状结构，并逐渐增大，当脉冲长度为200ms时，声辐射力太大破坏空化泡群枝状结构。图7给出脉冲长度为20ms时，阵列换能器位置从-5到6mm的各个空化泡群断面分布。图7给出脉冲长度为20ms时，阵列换能器位置从-5到6mm的空间序列空化分布。图8给出了脉冲长度为20ms时，三维空化图像(a)和三维空化微泡密度图像(b)。
实施例2.脉动流条件下，高强度聚焦超声脉冲长度20us，电功率100W，脉动频率80次/秒，流速为5cm/s，管道直径为4mm，对该实验条件下产生的空化泡群进行空化以及空化密度定量三维成像：
(1)制备琼脂明胶仿血管组织体模，将琼脂(3％)、明胶(12％)和除气水(85％)混合加热搅拌至完全溶解，放入除气箱除气30分钟，倒入模具待冷却凝固后慢慢抽出管子，形成管道，管道两端与乳胶管相连，并连接到脉动泵，脉动泵将液体从储存池经管道流入废液池，在管道内形成脉动流。
(2)采用图4的实验装置，水槽中注入自来水，高强度聚焦超声换能器固定于水槽侧壁，采用Sonix-Touch全数字化超声设备发射宽波束探测空化微泡，阵列换能器垂直于高强度聚焦超声传播方向放置，琼脂明胶仿血管组织体模与三维移动装置相连，通过三维移动将高强度聚焦超声聚焦在管道内部，在管道内产生空化。
(3)脉动频率设置为80次/秒，流速为5cm/s，储存液为生理盐水造影微泡溶液，将Sonix-Touch全数字化超声设备设置为脉冲多普勒模式，记录该脉动参数下的多普勒谱，分析其脉动周期和规律，以此编辑波形发生器AWG420的波形，其CH1、CH2和Marker通道分别输入给功率放大器、Sonix-Touch和脉动泵，实现三者之间的同步。
(4)固定阵列换能器位置，采用阵列Plane-by-plane宽波束检测空化，改变脉动泵与功率放大器之间的时间延迟，分别设置为0-12ms，单元时间为1ms，每次停止2s使得液体空化核分布恢复初始状态，得脉动周期内各时刻的序列二维空化分布图像原始射频数据。
(5)通过三维机械扫描装置Multiscan 5800沿垂直于阵列放置方向移动，移动单元距离设置为1mm，得到沿着得到不同单元位置的序列二维空化原始射频数。
(6)采用最小方差自适应波束合成算法对脉动周期内各时刻的序列二维空化分布图像原始射频数据和序列二维空化原始射频数据进行波束合成，得到脉动周期内各时刻和不同单元位置的序列二维空化射频数据。
(7)针对序列二维空化射频数据，一方面采用包络检波、对数压缩以及坐标变换得到脉动周期内各时刻和不同单元位置的序列二维空化分布图像，一方面采用Nakagami参量成像算法得到脉动周期内各时刻和不同单元位置的序列二维空化密度分布图像；
(8)采用三维重建算法，分别对序列二维空化分布图像和空化密度分布图像进行三维重建，得到高强度聚焦超声脉冲长度20us，电功率100W时脉动流下管道内三维空化分布图像和三维空化密度图像。
图9给实施例2条件下三维空化图像(a)和三维空化微泡密度图像(b)，可见空化微泡充满了仿血管组织体模管道，并沿着管道两端延伸。
使用宽波束检测空化后，阵列换能器移动一个单元位置，待空化核分布恢复，再在同样的空化能量激励下，再次宽波束检测空化，得到不同单元位置的空间系列二维空化原始射频数据，通过改变空化源能量作用时间、能量源装置激励与阵列换能器发射宽波束之间的时间延迟以及脉动泵与能量源装置激励之间的时间延迟，可得到时间序列二维空化原始射频数据，然后结合宽波束最小方差自适应波束合成、Nakagami参量成像、三维重建算法得到微秒分辨空化三维时空分布图像和空化微泡密度的定量图像，该方法结合阵列Plane-by-plane宽波束检测空化，宽波束最小方差自适应波束合成以及Nakagami参量算法成像，不仅具有空间高分辨、高信噪比的成像特点，而且可实现稳态自由场和脉动流条件下微秒分辨空化三维时空分布成像和空化微泡密度的定量成像，有潜力发展成为类似于声场测量的空化成像标准方法。

Claims (7)

  1. 微秒分辨空化时空分布的三维空化定量成像方法，其特征在于：包括以下步骤：采用宽波束检测空化得到二维空化原始射频数据，每次宽波束检测空化后，将用于进行宽波束检测的阵列换能器沿垂直于该阵列换能器放置方向移动一个单元位置，待空化核分布恢复到初始状态，再在同样的空化能量激励下，再次采用宽波束检测对应时间的空化，以此逐步得到对应该阵列换能器不同放置位置的一系列二维空化原始射频数据，然后再结合宽波束最小方差自适应波束合成、Nakagami参量成像和三维重建算法，得到三维空化图像和空化微泡密度的定量三维图像。
  2. 根据权利要求1所述微秒分辨空化时空分布的三维空化定量成像方法，其特征在于：空化的产生装置包括能量源装置和同步信号发生器；空化的检测装置包括可编程的发射宽波束的阵列换能器和并行通道数据采集及存储单元；同步信号发生器产生同步信号分别控制所述能量源装置和阵列换能器，能量源装置产生能量激励空化的产生，且能量源装置产生的能量连续可调，所述阵列换能器发射宽波束对空化进行检测，得到的空化回波信号由并行通道数据采集及存储单元采集存储。
  3. 根据权利要求1所述微秒分辨空化时空分布的三维空化定量成像方法，其特征在于：采用宽波束最小方差自适应波束合成算法对所述二维空化原始射频数据进行波束合成，得到沿垂直于所述阵列换能器放置方向不同放置位置的一系列二维空化射频数据，根据所述二维空化射频数据获得一系列二维空化分布图像，利用三维重建算法将所述二维空化分布图像进行三维显示，得到三维空化图像。
  4. 根据权利要求1所述微秒分辨空化时空分布的三维空化定量成像方法，其特征在于：采用宽波束最小方差自适应波束合成算法对所述二维空化原始射频数据进行波束合成，得到沿垂直于所述阵列换能器放置方向不同放置位置的一系列二维空化射频数据，对所述二维空化射频数据进行Nakagami参量提取后得到一系列二维空化密度定量图像，利用三维重建算法将所述二维空化密度定量图像进行三维显示，得到空化微泡密度的定量三维图像。
  5. 根据权利要求4所述微秒分辨空化时空分布的三维空化定量成像方法，其特征在于：所述Nakagami参量提取包括以下步骤：
    1)将二维空化射频数据rf进行去噪处理：
    a)选取一定尺寸大小的背景信号区域，计算该区域平均能量P；
    b)对rf分别叠加平均能量为P的随机高斯白噪声n1和n2，得到S1和S2；
    c)计算S1和S2的相关系数，并给定阈值Th，然后对相关系数进行阈值化处理后与rf加权，得到rfdenoise：
    Figure PCTCN2015071540-appb-100001
    ，其中corrcoef(S1,S2)为S1和S2的相关系数；
    2)将rfdenoise进行Hilbert解调，得到包络信号R；
    3)计算Nakagami参量：
    Figure PCTCN2015071540-appb-100002
    其中E(·)代表统计均值。
  6. 根据权利要求1所述微秒分辨空化时空分布的三维空化定量成像方法，其特征在于：在所述阵列换能器的同一个放置位置，根据脉动流多普勒谱记录的脉动周期和脉动规律，通过改变空化能量激励与脉动起始之间的时间延迟，得到脉动周期内不同时刻对应该放置位置的二维空化原始射频数据，在所述阵列换能器移动一个单元位置后，得到脉动周期内不同时刻对应新的放置位置的二维空化原始射频数据，最终得到脉动周期内任一时刻对应所述陈列换能器不同放置位置的一系列二维空化原始射频数据。
  7. 根据权利要求1所述微秒分辨空化时空分布的三维空化定量成像方法，其特征在于：改变空化能量激励时间，得到随时间演化的三维空化图像和空化微泡密度的定量三维图像；改变空化能量大小，得到随能量变化的三维空化图像和空化微泡密度的定量三维图像；改变空化能量激励与所述阵列换能器发射宽波束之间的时间延迟，得到随时间消散的三维空化图像和空化微泡密度的定量三维图像；改变空化能量作用的媒介，得到不同媒介下三维空化图像和空化微泡密度的定量三维图像。
PCT/CN2015/071540 2014-12-27 2015-01-26 微秒分辨空化时空分布的三维空化定量成像方法 WO2016101382A1 (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
US15/115,682 US10332250B2 (en) 2014-12-27 2015-01-26 Three-dimensional cavitation quantitative imaging method for microsecond-resolution cavitation spatial-temporal distribution

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
CN201410834392.1 2014-12-27
CN201410834392.1A CN104535645B (zh) 2014-12-27 2014-12-27 微秒分辨空化时空分布的三维空化定量成像方法

Publications (1)

Publication Number Publication Date
WO2016101382A1 true WO2016101382A1 (zh) 2016-06-30

Family

ID=52851215

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/CN2015/071540 WO2016101382A1 (zh) 2014-12-27 2015-01-26 微秒分辨空化时空分布的三维空化定量成像方法

Country Status (3)

Country Link
US (1) US10332250B2 (zh)
CN (1) CN104535645B (zh)
WO (1) WO2016101382A1 (zh)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116115260A (zh) * 2023-01-03 2023-05-16 西安交通大学 高分辨高对比快速计算的发收时序同步超声被动空化成像方法及系统

Families Citing this family (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103813752B (zh) * 2012-01-27 2017-11-10 东芝医疗系统株式会社 医用图像处理装置
CN104887266B (zh) * 2015-05-29 2017-04-26 西安交通大学 基于面阵的小区域三维被动空化成像及三维复合成像方法
CN107260217B (zh) * 2017-07-17 2018-07-17 西安交通大学 用于脑部聚焦超声空化实时监控的三维无源成像方法及系统
CN109431536B (zh) * 2018-09-17 2019-08-23 西安交通大学 一种聚焦超声空化的实时高分辨时空分布成像方法与系统
CN109523537B (zh) * 2018-11-20 2022-12-09 北京理工大学 一种基于图像处理的非定常空化流动结构精细化分析方法
CN109513123B (zh) * 2018-12-28 2020-05-19 西安交通大学 一种基于半球阵的高分辨三维被动空化成像方法
CN110206743B (zh) * 2019-05-28 2020-05-08 浙江大学 一种基于噪声纹理和气泡形态对照的轴流泵空化特征提取方法
CN111220700A (zh) * 2019-12-09 2020-06-02 中北大学 超声空化泡运动矢量估计方法
CN111000539A (zh) * 2019-12-30 2020-04-14 北京丰顺科技发展有限公司 脉象采集、诊脉装置以及脉象采集数据处理方法
CN113484470B (zh) * 2021-07-01 2023-07-04 浙江始祖鸟环境工程有限公司 一种高灵敏度的硫化氢气体检测器

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20050254345A1 (en) * 2000-07-10 2005-11-17 Sez America, Inc. Method and device for measuring cavitation
WO2013081850A1 (en) * 2011-11-29 2013-06-06 Nalco Company Fouling reduction device and method
CN103235041A (zh) * 2013-04-26 2013-08-07 西安交通大学 基于超声主动空化成像的空化起始阈值分布重建方法
CN103860201A (zh) * 2014-03-06 2014-06-18 西安交通大学 一种基于宽波束造影成像及提取灌注时间强度曲线的方法

Family Cites Families (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN1847824B (zh) * 2006-02-27 2011-01-05 西安交通大学 基于超高速摄影技术的超声场中微泡行为分析系统和方法
US8343050B2 (en) * 2009-05-04 2013-01-01 Siemens Medical Solutions Usa, Inc. Feedback in medical ultrasound imaging for high intensity focused ultrasound
RU2014116967A (ru) * 2011-09-27 2015-11-10 Конинклейке Филипс Н.В. Высокоинтенсивный фокусированный ультразвук, усиленный с помощью кавитации
JP6450752B2 (ja) * 2013-06-28 2019-01-09 コーニンクレッカ フィリップス エヌ ヴェＫｏｎｉｎｋｌｉｊｋｅ Ｐｈｉｌｉｐｓ Ｎ．Ｖ． 画像誘導超音波血栓溶解のためのトランスデューサ配置及び位置合わせ

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20050254345A1 (en) * 2000-07-10 2005-11-17 Sez America, Inc. Method and device for measuring cavitation
WO2013081850A1 (en) * 2011-11-29 2013-06-06 Nalco Company Fouling reduction device and method
CN103235041A (zh) * 2013-04-26 2013-08-07 西安交通大学 基于超声主动空化成像的空化起始阈值分布重建方法
CN103860201A (zh) * 2014-03-06 2014-06-18 西安交通大学 一种基于宽波束造影成像及提取灌注时间强度曲线的方法

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116115260A (zh) * 2023-01-03 2023-05-16 西安交通大学 高分辨高对比快速计算的发收时序同步超声被动空化成像方法及系统
CN116115260B (zh) * 2023-01-03 2024-05-24 西安交通大学 高分辨高对比快速计算的发收时序同步超声被动空化成像方法及系统

Also Published As

Publication number Publication date
US10332250B2 (en) 2019-06-25
US20170011508A1 (en) 2017-01-12
CN104535645B (zh) 2016-06-29
CN104535645A (zh) 2015-04-22

Similar Documents

Publication Publication Date Title
WO2016101382A1 (zh) 微秒分辨空化时空分布的三维空化定量成像方法
CN110248606B (zh) 气穴定位
OˈReilly et al. A super‐resolution ultrasound method for brain vascular mapping
Salgaonkar et al. Passive cavitation imaging with ultrasound arrays
Xu et al. Effects of acoustic parameters on bubble cloud dynamics in ultrasound tissue erosion (histotripsy)
Hall et al. A real-time measure of cavitation induced tissue disruption by ultrasound imaging backscatter reduction
Song et al. Coded excitation plane wave imaging for shear wave motion detection
KR20180013956A (ko) 단일 추적 위치 전단파 탄성 이미징을 위한 방법, 시스템 및 컴퓨터 프로그램 제품
US11571230B2 (en) Ultrasonic processing apparatus comprising means for imaging cavitation bubbles
Xu et al. Optical and acoustic monitoring of bubble cloud dynamics at a tissue-fluid interface in ultrasound tissue erosion
KR101610874B1 (ko) 공간 일관성 기초 초음파 신호 처리 모듈 및 그에 의한 초음파 신호 처리 방법
Ding et al. Ultrasound line-by-line scanning method of spatial–temporal active cavitation mapping for high-intensity focused ultrasound
Hu et al. Spatial-temporal ultrasound imaging of residual cavitation bubbles around a fluid–tissue interface in histotripsy
US20090187107A1 (en) Ultrasonic diagnostic apparatus and control program thereof
Miller et al. Bubble-induced color Doppler feedback for histotripsy tissue fractionation
Li et al. Tri-modality cavitation mapping in shock wave lithotripsy
Osipov et al. 3D ultrasound: Current state, emerging trends and technologies
Ding et al. Spatial–temporal three-dimensional ultrasound plane-by-plane active cavitation mapping for high-intensity focused ultrasound in free field and pulsatile flow
Park et al. Photoacoustic imaging using array transducer
Filoux et al. High-frequency annular array with coaxial illumination for dual-modality ultrasonic and photoacoustic imaging
Collin et al. Real-time three-dimensional passive cavitation detection for clinical high intensity focussed ultrasound systems
O'Reilly et al. Investigating a method for non-invasive ultrasound aberration correction through the skull bone
Kim et al. Comparison study of passive acoustic mapping and high-speed photography for monitoring in situ cavitation bubbles
CN111220700A (zh) 超声空化泡运动矢量估计方法
US9247921B2 (en) Systems and methods of high frame rate streaming for treatment monitoring

Legal Events

Date Code Title Description
WWE Wipo information: entry into national phase

Ref document number: 15115682

Country of ref document: US

121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 15871461

Country of ref document: EP

Kind code of ref document: A1

NENP Non-entry into the national phase

Ref country code: DE

122 Ep: pct application non-entry in european phase

Ref document number: 15871461

Country of ref document: EP

Kind code of ref document: A1