CN110772281B - Ultrasonic CT imaging system based on improved ray tracing method - Google Patents
Ultrasonic CT imaging system based on improved ray tracing method Download PDFInfo
- Publication number
- CN110772281B CN110772281B CN201911012719.6A CN201911012719A CN110772281B CN 110772281 B CN110772281 B CN 110772281B CN 201911012719 A CN201911012719 A CN 201911012719A CN 110772281 B CN110772281 B CN 110772281B
- Authority
- CN
- China
- Prior art keywords
- sound
- ultrasonic
- imaging
- target
- organ
- Prior art date
- Legal status (The legal status 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 status listed.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 70
- 238000013170 computed tomography imaging Methods 0.000 title claims abstract description 22
- 238000003384 imaging method Methods 0.000 claims abstract description 50
- 210000000056 organ Anatomy 0.000 claims abstract description 40
- 238000012545 processing Methods 0.000 claims abstract description 10
- 238000004590 computer program Methods 0.000 claims abstract description 9
- 238000002604 ultrasonography Methods 0.000 claims description 31
- 238000004422 calculation algorithm Methods 0.000 claims description 15
- 230000008569 process Effects 0.000 claims description 8
- 230000015572 biosynthetic process Effects 0.000 claims description 6
- 238000003786 synthesis reaction Methods 0.000 claims description 6
- 238000001914 filtration Methods 0.000 claims description 5
- 239000000523 sample Substances 0.000 claims description 4
- 238000005070 sampling Methods 0.000 claims description 2
- 210000000481 breast Anatomy 0.000 description 35
- 238000010586 diagram Methods 0.000 description 6
- 238000004364 calculation method Methods 0.000 description 5
- 239000011159 matrix material Substances 0.000 description 5
- 230000005236 sound signal Effects 0.000 description 4
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 description 4
- 230000005540 biological transmission Effects 0.000 description 3
- 230000001276 controlling effect Effects 0.000 description 3
- 238000004088 simulation Methods 0.000 description 3
- 230000003111 delayed effect Effects 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 230000003028 elevating effect Effects 0.000 description 2
- 238000011156 evaluation Methods 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 230000001902 propagating effect Effects 0.000 description 2
- 230000001105 regulatory effect Effects 0.000 description 2
- 238000012285 ultrasound imaging Methods 0.000 description 2
- 206010006187 Breast cancer Diseases 0.000 description 1
- 208000026310 Breast neoplasm Diseases 0.000 description 1
- 238000004891 communication Methods 0.000 description 1
- 239000007822 coupling agent Substances 0.000 description 1
- 238000000605 extraction Methods 0.000 description 1
- 238000011478 gradient descent method Methods 0.000 description 1
- 238000007689 inspection Methods 0.000 description 1
- 238000012804 iterative process Methods 0.000 description 1
- 238000013178 mathematical model Methods 0.000 description 1
- 239000002245 particle Substances 0.000 description 1
- 230000005855 radiation Effects 0.000 description 1
- 238000012216 screening Methods 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
- 230000007704 transition Effects 0.000 description 1
- 238000012795 verification Methods 0.000 description 1
Images
Classifications
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B8/00—Diagnosis using ultrasonic, sonic or infrasonic waves
- A61B8/08—Detecting organic movements or changes, e.g. tumours, cysts, swellings
- A61B8/0825—Detecting organic movements or changes, e.g. tumours, cysts, swellings for diagnosis of the breast, e.g. mammography
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B8/00—Diagnosis using ultrasonic, sonic or infrasonic waves
- A61B8/13—Tomography
- A61B8/15—Transmission-tomography
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B8/00—Diagnosis using ultrasonic, sonic or infrasonic waves
- A61B8/52—Devices using data or image processing specially adapted for diagnosis using ultrasonic, sonic or infrasonic waves
- A61B8/5207—Devices using data or image processing specially adapted for diagnosis using ultrasonic, sonic or infrasonic waves involving processing of raw data to produce diagnostic data, e.g. for generating an image
Landscapes
- Health & Medical Sciences (AREA)
- Life Sciences & Earth Sciences (AREA)
- Engineering & Computer Science (AREA)
- Biomedical Technology (AREA)
- Medical Informatics (AREA)
- Pathology (AREA)
- Radiology & Medical Imaging (AREA)
- Biophysics (AREA)
- Physics & Mathematics (AREA)
- Heart & Thoracic Surgery (AREA)
- Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
- Molecular Biology (AREA)
- Surgery (AREA)
- Animal Behavior & Ethology (AREA)
- General Health & Medical Sciences (AREA)
- Public Health (AREA)
- Veterinary Medicine (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Ultra Sonic Daignosis Equipment (AREA)
Abstract
The invention relates to an ultrasonic CT imaging system based on an improved ray tracing method, which comprises a full-surrounding ultrasonic transducer, a signal acquisition module and a data processing and imaging module. The processor in the data processing and imaging module, when executing the computer program stored in the memory, performs the steps of: scanning human organs based on a fully-enclosed ultrasonic transducer, and measuring the transit time of sound waves; reconstructing initial sound velocity distribution of a target; b ultrasonic imaging is carried out based on the sound velocity distribution of the target; extracting organ tissue boundary information from the B ultrasonic image; iteratively updating the propagation path of the sound wave; the sound speed distribution of the target is iteratively updated. Compared with the traditional imaging scheme, the invention reduces the calculated amount, improves the reconstruction speed and improves the precision of the B-mode ultrasonic image and the edge definition in the ultrasonic CT imaging result while ensuring the ultrasonic CT imaging precision.
Description
Technical Field
The invention relates to an ultrasonic CT imaging system, in particular to an improved ray tracing method and a system for carrying out ultrasonic CT imaging by utilizing breast tissue boundary information extracted from B-mode ultrasound based on an algorithm.
Background
The basic principle of the ultrasonic CT as a new ultrasonic tomographic image reconstruction method is to scan a breast by using ultrasonic waves and extract effective information such as transit time from an ultrasonic signal to reconstruct acoustic properties of breast tissue such as sound velocity distribution. The method has the characteristics of high resolution, no radiation and the like, and has wide application prospects in the fields of nondestructive inspection, early breast cancer screening and the like.
A conventional ultrasound CT imaging procedure can be divided into three steps: scanning the breast with a fully-enclosed probe and extracting the transit time from the acoustic signal; reconstructing a propagation path of sound in breast tissue; based on the transit time and the acoustic propagation path, the sound velocity distribution in the imaging region is solved.
In the above steps, it is a difficult point to reconstruct an accurate sound propagation path. The reconstruction methods of ultrasound CT can be divided into two categories according to the method of reconstructing the acoustic path. In the first method, the ultrasonic wave propagates along a straight line, and the propagation path is a straight line from the transmitting array element to the receiving array element. The method is simple in calculation, but phenomena such as refraction and reflection of sound in mammary tissue are not considered, so that the reconstruction result based on the method has the problems of boundary blurring, shape distortion and the like. The second method models the propagation process of sound more accurately, so that the ultrasonic wave propagates along a curve. In the method, firstly, the transit time of sound waves reaching each point in an imaging area is calculated based on the Huygens principle, and then a route with the shortest transit time is found as a sound path by using a gradient descent method, and the method is called a ray tracing method. Compared with the first method, the ray tracing method considers the refraction and diffraction phenomena of sound in the transmission process, and the mathematical model of the ray tracing method is closer to the real model, so that the reconstruction result is more accurate; the disadvantage is that the acoustic propagation path needs to be optimized step by step in an iterative process, and the required amount of computation is extremely large, which reduces the practicality of the method.
Compared with the traditional ultrasonic imaging (namely B-ultrasonic imaging) mode, the reconstruction result of the ultrasonic CT can intuitively reflect the sound velocity distribution of the mammary tissue. However, the imaging method is not sensitive to the boundary due to the reconstruction principle, and the reconstruction result often has the problems of boundary blurring and large calculation amount. On the contrary, the B-mode ultrasonic imaging utilizes the echo to detect the boundary of the mammary tissue, the reconstruction time is short, the boundary is obvious, but the speckle noise in the reconstruction result is serious. In addition, there are several assumptions in B-ultrasonic imaging, such as tissue sound velocity of 1540m/s, ultrasonic wave propagating along a straight line, etc., and these assumptions are not strictly true in the actual reconstruction process, so there is often serious shape distortion in the reconstruction result.
The invention provides a technical scheme which is based on an improved ray tracing method and utilizes mammary tissue boundary information extracted from B ultrasonic to calculate the propagation path of ultrasonic waves in a breast. The ultrasonic propagation path obtained by the scheme is close to the result of the traditional ray tracing method, but the calculation amount is reduced by more than 90%. The ultrasonic CT is reconstructed based on the algorithm, so that the reconstruction time can be reduced by more than one order of magnitude while the high precision of the reconstruction result is ensured; in addition, the invention can simultaneously improve the precision of the B ultrasonic image and the edge definition in the ultrasonic CT imaging result.
Disclosure of Invention
The invention aims to improve the traditional ray tracing method and use the breast tissue boundary information extracted from the B-ultrasonic image for ultrasonic CT imaging based on the algorithm. Compared with the ultrasonic CT imaging based on the traditional ray tracing method, the method has the advantages that the calculation amount is reduced, the reconstruction speed is increased, and the precision of the B ultrasonic image and the edge definition in the ultrasonic CT imaging result are improved while the ultrasonic CT imaging precision is ensured.
The first aspect of the technical solution of the present invention is an ultrasound CT imaging system, comprising: a fully-enclosed ultrasound transducer configured to surround the determined location of the human organ;
the signal acquisition module is used for sampling the acoustic signals received by the all-around ultrasonic transducer;
a data processing and imaging module coupled to the signal acquisition module, the data processing and imaging module including a memory and a processor, wherein the processor implements a method when executing a computer program stored in the memory, the method comprising:
scanning human organs based on a fully-enclosed ultrasonic transducer, and measuring the transit time of sound waves;
step two, reconstructing initial sound velocity distribution of the target;
step three, B ultrasonic imaging is carried out based on the sound velocity distribution of the target;
step four, extracting organ tissue boundary information from the B ultrasonic image;
step five, iteratively updating the propagation path of the sound wave;
and step six, iteratively updating the sound velocity distribution of the target.
In some embodiments, the first step comprises: the fully-enclosed ultrasonic transducer is caused to transmit and receive ultrasonic waves and the time of flight is extracted from the recorded ultrasonic signal based on a "threshold-phase" method.
In some embodiments, the second step comprises: the propagation path of the recording sound in the organ tissue is a straight line from the transmitting array element to the receiving array element, and is L(0)Represents; recording the slowness of sound as s; and solving formula L based on SART algorithm(0)s(0)Get an initial solution s of s ═ t(0)。
In some embodiments, the third step comprises: based on the full-enclosure ultrasonic transducer, digital beam synthesis is carried out in a mode of electronic focusing and delay summation, and B ultrasonic imaging is carried out on human organs; and in the process of delay summation, calculating the emission delay time of each emission array element based on the sound slowness distribution of the organ tissues reconstructed in the second step or the fifth step.
In some embodiments, the fourth step comprises: performing convolution and operation on each pixel point in the B ultrasonic imaging result by using a Sobel operator, and then selecting a threshold value to extract an edge; after the filtering operation is performed, the boundary of the organ tissue is generated.
In some embodiments, said step five comprises: modeling refraction and reflection phenomena occurring on organ tissue boundaries during acoustic propagation based on acoustic slowness s(r)Calculating the sound arrival time of the ultrasonic wave reaching each point of the organ tissue boundary, and iteratively updating the propagation path of the ultrasonic wave in the organ tissue to obtain L(r+1)。
In some embodiments, the sixth step comprises: iteratively updating the sound velocity distribution of the target, and solving based on SART algorithmL(r+1)s(r+1)Obtaining the updated sound slowness s as t(r+1)。
In some embodiments, the fully-enclosed ultrasonic transducer comprises a ring-shaped ultrasonic probe composed of a plurality of array elements. The system also includes a lifting device, wherein the lifting device includes: the supporting frame is used for supporting the whole fully-enclosed ultrasonic transducer and the attached scanning equipment; an inner cavity filled with an ultrasound couplant to receive a body organ; and the motor is used for driving the all-around ultrasonic transducer to perform linear motion along the normal direction of the reference plane of the human organ. The system also comprises a control module which controls the emission and the reception of sound waves of the fully-enclosed ultrasonic transducer through the switching of an electronic switch on one hand and adjusts the height of the transducer through controlling a lifting device on the other hand so as to realize the three-dimensional scanning of human organs
A second aspect of the present invention is a computer-readable storage medium having a computer program stored thereon. The computer program when executed by a processor implements a method comprising:
scanning human organs based on a fully-enclosed ultrasonic transducer, and measuring the transit time of sound waves;
step two, reconstructing initial sound velocity distribution of the target;
step three, B ultrasonic imaging is carried out based on the sound velocity distribution of the target;
step four, extracting organ tissue boundary information from the B ultrasonic image;
step five, iteratively updating the propagation path of the sound wave;
and step six, iteratively updating the sound velocity distribution of the target.
Compared with the prior art, the invention has the following advantages:
the traditional ray tracing method is improved, the propagation path and the sound slowness distribution of the sound wave are reconstructed on the basis of the mammary tissue boundary extracted from the B ultrasonic image, and compared with the traditional ray tracing method, the calculated amount is reduced by more than 90% while the ultrasonic CT imaging precision is ensured; meanwhile, the breast tissue boundary extracted from the B ultrasonic image is applied to the ultrasonic CT, so that the imaging result boundary of the scheme is clearer; in addition, electronic focusing and beam synthesis are carried out on the basis of the acoustic slowness distribution reconstructed by the ultrasonic CT, and the imaging result of the B ultrasonic is improved.
The imaging scheme provided by the invention has strong applicability and good imaging effect. Although the scheme is provided based on a fully-enclosed ultrasonic transducer, the scheme is not necessary, and as long as the ultrasonic transducer with an accurately known position is used for acquiring enough transmission type projection data for a human breast, for example, a plurality of linear ultrasonic transducers are used for rotating around the human breast and acquiring ultrasonic signals, ultrasonic CT imaging can be performed according to the reconstruction method provided by the invention; the scheme improves the imaging results of ultrasonic CT and B-ultrasonic simultaneously, wherein compared with the traditional scheme, the imaging result of ultrasonic CT reduces the calculated amount by more than 90 percent on the basis of ensuring that the reconstruction precision is not reduced.
Drawings
FIG. 1 is a block diagram of the system architecture of the present invention.
FIG. 2 is a flow chart of a method of implementing the system of the present invention.
Fig. 3 is a schematic diagram of a fully-enclosed ultrasound transducer in the system of the present invention for ultrasound CT imaging.
Fig. 4 is a schematic diagram of the extraction of transit time from an acoustic signal in accordance with the present invention.
Fig. 5 is a schematic diagram of the sound wave propagating along a straight line or a broken line in the present invention.
Fig. 6 is a simulated digital standard phantom used in the present invention to reconstruct the sound velocity distribution of breast tissue.
Fig. 7 is a schematic diagram of the initial value of the sound velocity distribution of the breast tissue obtained in the present invention.
Fig. 8 is a schematic diagram of the sound velocity distribution of the breast tissue obtained by the final reconstruction in the present invention.
Detailed Description
The following examples are provided to illustrate specific embodiments of the present invention. The ultrasound CT imaging system of the present invention can be applied to imaging of organs or tissues of the human body, and the imaging scheme is mainly described herein by taking the breast as an example, but not as a limitation of the imaging object.
The hardware architecture of the ultrasonic CT imaging system can comprise: elevating gear, surround ultrasonic transducer entirely, signal acquisition module, control module, data processing and imaging module, as shown in figure 1. Wherein, elevating gear includes support frame, interior cavity and motor. The fully-enclosed ultrasound transducer position may comprise a ring-shaped ultrasound probe consisting of a plurality of array elements, and its spatial position is predetermined. The plurality of array elements may be arranged in a circular, fully-enclosed manner, as shown in fig. 3 and 5. The support frame is used for supporting the whole scanning device (comprising a full-enclosure ultrasonic transducer), the inner cavity is used for receiving human breasts and containing an ultrasonic coupling agent (such as water), and the motor is connected with the full-enclosure ultrasonic transducer and used for controlling the lifting motion of the transducer in a vertical plane. The full-surrounding ultrasonic transducer surrounds the inner cavity of the lifting device and is used for transmitting and receiving ultrasonic signals. The data acquisition module samples the acoustic signals received by the all-around ultrasonic transducer and outputs the result to the data processing and imaging module. The control module controls the sound wave emission and the sound wave reception of the fully-enclosed ultrasonic transducer through the switching of the electronic switch on one hand, and adjusts the height of the transducer through controlling the lifting device on the other hand, so that the three-dimensional scanning of the human breast is realized. The data processing and imaging module is responsible for analyzing the acquired ultrasonic signals and finishing the reconstruction of the ultrasonic CT image and the B ultrasonic image.
Referring to fig. 2, the method for implementing ultrasound CT imaging by the data processing and imaging module mainly includes the following steps: firstly, scanning a human breast based on a fully-enclosed ultrasonic transducer, and measuring the transit time of sound waves; secondly, reconstructing the initial sound velocity distribution of the target; thirdly, B-ultrasonic imaging is carried out based on the sound velocity distribution of the target; fourthly, extracting the boundary information of the mammary tissue from the B-ultrasonic image; fifthly, iteratively updating the propagation path of the sound wave; and finally, iteratively updating the sound velocity distribution of the target. The following describes each step and operation process in detail.
The method comprises the following steps: scanning human breast based on full-enclosure ultrasonic transducer and measuring transit time of sound wave
In this step, a fully-enclosed ultrasound transducer is used to encircle the human breast, transmit and receive ultrasound waves, and the transit time is extracted from the recorded ultrasound signals based on a 'threshold-phase' method.
The method comprises the following specific steps:
(1.1) placing a human breast in the center of a fully-enclosed ultrasonic transducer with an array element number of N, and immersing the transducer and the breast in water to enable each array element of the transducer to sequentially emit sound waves; when each array element transmits sound wave, all array elements receive and record sound signals, and therefore total NxN sound signals are recorded;
(1.2) let the received acoustic signal be denoted by p ═ f (τ), where p denotes the sound pressure and τ denotes the propagation time of the sound; the corresponding transit time t of a single ultrasound signal is then given by:
where Δ p is a threshold value, which is slightly higher than the noise maximum in the recorded ultrasound signal, τ-And τ+Respectively representing a left neighborhood and a right neighborhood of tau; the meaning of the formula is: after the received ultrasonic signal is higher than the threshold value for the first time, taking a first maximum value point as a characteristic point P, wherein the time corresponding to P is the sound arrival time t of the ultrasonic signal; extracting the transit time from each ultrasonic signal to form an N2Vector t of x 1:
wherein, thRepresenting the time of flight extracted from the acoustic signal corresponding to the h-th transmit-receive array element pair.
Step two: reconstructing initial sound velocity distribution of target
The path of sound propagation in the mammary tissue is marked as a straight line from the transmitting array element to the receiving array element, denoted by L(0)Represents; the note slowness (i.e., the reciprocal of the speed of sound) is s; solving formula L based on SART algorithm(0)s(0)Get an initial solution s of s ═ t(0)。
The method comprises the following specific steps:
(2.1) to getThe image area is divided into M × M grids, the transmitting array element and the receiving array element are connected by a straight line, and the length of the straight line in each grid is obtained, because N is total2A transmit-receive array element pair, so that the result is N2×M2Matrix of dimensions, denoted
Where h denotes the h-th transmit-receive array element pair, k denotes the k-th trellis,
represents the length of the h-th acoustic path in the k-th grid;
(2.2) reconstructing the sound speed distribution in the imaging region: firstly, the speed of sound sinitIs initialized to the sound velocity of water, and then L is solved based on SART algorithm(0)s(0)＝t：
In the formula, λ is a parameter for regulating convergence speed; the solution of the above formula is in the form:
wherein,
indicating the initial value of the acoustic slowness in the kth grid.
Step three: b-mode ultrasound imaging based on target sound velocity distribution
In the step, digital beam synthesis is carried out by means of electronic focusing and delay summation based on a full-enclosure ultrasonic transducer, and B-ultrasonic imaging is carried out on human breasts; in the process of time delay summation, the emission delay time of each array element is calculated based on the sound slowness distribution of the mammary tissue reconstructed in the second step or the fifth step:
wherein Tr represents an array element for transmitting ultrasonic waves, Fo represents a focal point for focusing electrons,
is the delayed emission time of Tr when the focal point Fo is located, dTr,FoIs the length of the straight line from Tr to Fo in each grid, s(r)Represents the acoustic slowness reconstructed in step two or step five: when r is 0, s(r)Representing the sound slowness obtained by reconstruction in the step two; when r ≠ 0, s(r)And (4) representing the sound slowness obtained by reconstruction in the step five, wherein r represents the iteration number, and 1 is added to the value of r every time of iteration in the step five.
Step four: extracting mammary tissue boundary information from B-ultrasonic image
In the step, a Sobel operator is used for performing convolution and operation on each pixel point in a B ultrasonic imaging result, then a proper threshold value is selected to extract an edge, and after simple filtering operation is performed, the boundary of the mammary tissue can be obtained.
Step five: iteratively updating the propagation path of an acoustic wave
Modeling refraction and reflection phenomena occurring on the boundary of mammary tissue during sound propagation based on the acoustic slowness s(r)Calculating the sound arrival time of the ultrasonic wave reaching each point of the breast tissue boundary, and iteratively updating the propagation path of the ultrasonic wave in the breast tissue to obtain L(r+1)。
For each array element Tr for transmitting ultrasonic wave, the array element for receiving ultrasonic wave is Tr', and the following steps are implemented:
(5.1) for each point P on the boundary of the mammary tissueACalculating the arrival P of the ultrasonic wave from TrAThe initial value of the sound arrival time is specifically calculated in the following way: if Tr and PAThe connecting line between them passes through the boundary of mammary tissue
Otherwise
Wherein,
indicating the arrival of ultrasonic waves from Tr to PAThe sound of (a) is given for a time,
represents the length of the straight line from Tr to Fo within each grid;
(5.2) iteratively updating the arrival of Tr at each point P on the breast tissue boundaryASound arrival time of
In the formula,
represents the updated sound arrival time, PBRepresenting the division of P on the boundary of mammary tissueAAt any point of the exterior of the body,
indicates the ultrasonic wave from PBTo PAThe sound arrival time of (2) is calculated by the following method: if PBAnd PAThe connecting line between them passes through the boundary of mammary tissue
Otherwise
When for each point P on the upper boundary of the breast tissueAAll are provided with
Then, the iteration is stopped, wherein, the delta mu is a preset small parameter;
(5.3) calculating the minimum acoustic propagation time μ from Tr to TrTr,Tr′：
And the acoustic path corresponding to the acoustic propagation time is the updated acoustic propagation path.
Iteratively updating the sound propagation path between each transmitting-receiving array element pair, and recording the updated sound wave propagation path matrix as L(r+1)。
Step six: iteratively updating a sound velocity distribution of a target
Solving for L based on SART algorithm(r+1)s(r+1)Obtaining the updated sound slowness s as t(r+1)：
Where Δ s is a predetermined small parameter, if s(r+1)Satisfy the iteration stop condition | s(r+1)-s(r)|<Δ s, then s(r +1)If not, the final output sound slowness distribution is input into the third step as prior information, B-ultrasonic imaging is carried out again until the obtained s(r+1)The condition is satisfied.
In the following, the technical scheme of the ultrasonic CT imaging system of the present invention is subjected to simulation verification and evaluation by combining with a specific calculation example and based on an MATLAB simulation platform. The digital phantom designed and used in the examples is shown in fig. 6.
Executing the step one: scanning human breast based on full-enclosure ultrasonic transducer and measuring transit time of sound wave
Scanning the body model by using a full-surrounding annular ultrasonic transducer, recording an ultrasonic signal, and extracting the transit time based on a threshold-phase method.
The method comprises the following specific steps:
(Z1.1) scanning the phantom with a fully-enclosed annular ultrasound transducer with an array element number of 256 (see fig. 3 or 5): when each array element of the transducer emits sound waves in sequence, all the array elements receive and record sound signals when each array element emits the sound waves, and therefore 65536 sound signals are recorded in total;
(Z1.2) the waveform of the received acoustic signal p ═ f (τ) is shown in fig. 4, where p represents the sound pressure, and in this embodiment p is normalized, i.e., the maximum value of the absolute value of the signal amplitude is 1, and τ represents the propagation time of the sound; calculating the corresponding transition time t of the ultrasonic signal according to the following formula:
where Δ p is a threshold value, which is slightly higher than the maximum noise value in the recorded ultrasound signal, in this example, 0.07, τ-And τ+Representing the left and right neighbourhoods of τ, respectively. Calculating by the above formula to obtain the corresponding transit time of the ultrasonic signal to be 1.825 microseconds; the transit time is extracted from each ultrasound signal to form a 65536 x 1 vector t.
And (5) executing the step two: reconstructing initial sound velocity distribution of target
The path of sound propagation in the mammary tissue is marked as a straight line from the transmitting array element to the receiving array element, denoted by L(0)Represents; the note slowness (i.e., the reciprocal of the speed of sound) is s; solving formula L based on SART algorithm(0)s(0)Get an initial solution s of s ═ t(0)。
The method comprises the following specific steps:
(Z2.1) the imaging region is divided into 256 × 256 grids, the transmission array elements and the reception array elements are connected by straight lines, and the length of the straight lines in each grid is determined, and the total number of the transmission-reception array element pairs is 65536 × 65536, so that the result is a 65536 × 65536 matrix, which is referred to as a "65536 × 65536" - "matrix"
Where h denotes the h-th transmit-receive array element pair, k denotes the k-th trellis,
represents the length of the h-th acoustic path in the k-th grid;
(Z2.2) reconstructing the sound speed distribution in the imaging region: firstly, the speed of sound sinitInitialized to the sound velocity of water 1515m/s, and then solved for L based on SART algorithm(0)s(0)＝t：
In the formula, λ is a parameter for regulating convergence speed, the solution of the above formula is a 65536 × 1 vector, each element represents the sound slowness value in 65536 grids, and then the velocity value distribution in the imaging region is reconstructed according to the sound slowness value, and the result is shown in fig. 7.
And step three is executed: b-mode ultrasound imaging based on target sound velocity distribution
In the step, digital beam synthesis is carried out by means of electronic focusing and delay summation based on a full-enclosure ultrasonic transducer, and B-ultrasonic imaging is carried out on human breasts; in the process of time delay summation, the emission delay time of each array element is calculated based on the sound slowness distribution of the mammary tissue reconstructed in the second step or the fifth step:
wherein Tr represents an array element for transmitting ultrasonic waves, Fo represents a focal point for focusing electrons,
is the delayed emission time of Tr when the focal point Fo is located, dTr,FoIs the length of the straight line from Tr to Fo in each grid, s(r)In step two or step fiveReconstructed acoustic slowness: when r is 0, s(r)Representing the sound slowness obtained by reconstruction in the step two; when r ≠ 0, s(r)And (4) representing the sound slowness obtained by reconstruction in the step five, wherein r represents the iteration number, and 1 is added to the value of r every time of iteration in the step five.
And step four is executed: extracting mammary tissue boundary information from B-ultrasonic image
In the step, a Sobel operator is used for performing convolution and operation on each pixel point in a B ultrasonic imaging result, then a proper threshold value is selected to extract an edge, and after simple filtering operation is performed, the boundary of the mammary tissue can be obtained.
And executing the step five: iteratively updating the propagation path of an acoustic wave
Modeling refraction and reflection phenomena occurring on the boundary of mammary tissue during sound propagation based on the acoustic slowness s(r)Calculating the sound arrival time of the ultrasonic wave reaching each point of the breast tissue boundary, and iteratively updating the propagation path of the ultrasonic wave in the breast tissue to obtain L(r+1)。
For each array element Tr for transmitting ultrasonic wave, the array element for receiving ultrasonic wave is Tr', and the following steps are implemented:
(Z5.1) for each point P on the boundary of the mammary tissueACalculating the arrival P of the ultrasonic wave from TrAThe initial value of the sound arrival time is specifically calculated in the following way: if Tr and PAThe connecting line between them passes through the boundary of mammary tissue
Otherwise
Wherein,
indicating the arrival of ultrasonic waves from Tr to PAThe sound of (a) is given for a time,
represents the length of the straight line from Tr to Fo within each grid;
(Z5.2) iteratively updating the arrival of Tr at each point P on the breast tissue boundaryASound arrival time of
In the formula,
represents the updated sound arrival time, PBRepresenting the division of P on the boundary of mammary tissueAAt any point of the exterior of the body,
indicates the ultrasonic wave from PBTo PAThe sound arrival time of (c); when for each point P on the upper boundary of the breast tissueAAll are provided with
When the iteration is stopped, wherein the Δ μ is 1/80 microseconds, which is a preset small parameter;
(Z5.3) calculating the minimum acoustic propagation time μ from Tr to TrTr,Tr′：
And the acoustic path corresponding to the acoustic propagation time is the updated acoustic propagation path.
Iteratively updating the sound propagation path between each transmitting-receiving array element pair, and recording the updated sound wave propagation path matrix as L(r+1)。
And a sixth step is executed: iteratively updating a sound velocity distribution of a target
Solving for L based on SART algorithm(r+1)s(r+1)Obtaining the updated sound slowness s as t(r+1)：
Existing Δ s is 2 × 10-6(m/s)-1Is a preset small parameter, if s(r+1)Satisfy the iteration stop condition | s(r+1)-s(r)|<Δ s, then s(r+1)If not, inputting the final output sound slowness distribution as prior information, executing the step three, and performing B-ultrasonic imaging again until the obtained s(r+1)Satisfy the condition, in this embodiment, s after one iteration(r+1)And (6) converging.
In this embodiment, the sound velocity distribution of the breast tissue reconstructed based on the improved ray tracing method is shown in fig. 8, and the reconstructed result is substantially equivalent to the result based on the conventional ray tracing method in quality: the root mean square error between the sound velocity distribution of the reconstruction result and the sound velocity distribution of the simulation phantom is used as an evaluation index, and the difference between the results of the two methods is not more than 2.5%; particularly, the boundary definition of a reconstruction result based on the algorithm provided by the invention is higher than that of the traditional ray tracing method by more than 50%; in addition, the sound velocity distribution of the breast tissue reconstructed based on the improved ray tracing method is reduced by more than 90% compared with the traditional ray tracing method, and the reconstruction time can be saved, in this embodiment, it takes 952 seconds to reconstruct the sound velocity distribution of the breast tissue based on the improved ray tracing method, while the traditional ray tracing method takes 12018 seconds.
It should be recognized that embodiments of the present invention can also be implemented or realized in computer hardware, a combination of hardware and software, or by computer instructions stored in a non-transitory computer readable memory. The methods may be implemented in a computer program using standard programming techniques, including a non-transitory computer readable storage medium configured with the computer program, where the storage medium so configured causes a computer to operate in a specific and predefined manner, according to the methods and figures described in the detailed description. Each program may be implemented in a high level procedural or object oriented programming language to communicate with a computer system. However, the program(s) can be implemented in assembly or machine language, if desired. In any case, the language may be a compiled or interpreted language. Furthermore, the program can be run on a programmed application specific integrated circuit for this purpose.
Further, the above-described methods may be implemented in any type of computing platform operatively connected to a suitable connection, including but not limited to a personal computer, mini computer, mainframe, workstation, networked or distributed computing environment, separate or integrated computer platform, or in communication with a charged particle tool or other imaging device, and the like. Aspects of the invention may be embodied in machine-readable code stored on a non-transitory storage medium or device, whether removable or integrated into a computing platform, such as a hard disk, optically read and/or write storage medium, RAM, ROM, or the like, such that it may be read by a programmable computer, which when read by the storage medium or device, is operative to configure and operate the computer to perform the procedures described herein. Further, the machine-readable code, or portions thereof, may be transmitted over a wired or wireless network. The invention described herein includes these and other different types of non-transitory computer-readable storage media when such media include instructions or programs that implement the steps described above in conjunction with a microprocessor or other data processor. The invention also includes the computer itself when programmed according to the methods and techniques described herein.
The above description is only a preferred embodiment of the present invention, and the present invention is not limited to the above embodiment, and any modifications, equivalent substitutions, improvements, etc. within the spirit and principle of the present invention should be included in the protection scope of the present invention as long as the technical effects of the present invention are achieved by the same means. The invention is capable of other modifications and variations in its technical solution and/or its implementation, within the scope of protection of the invention.
Claims (4)
1. An ultrasound CT imaging system comprising:
a fully-enclosed ultrasound transducer configured to surround the determined location of the human organ;
the signal acquisition module is used for sampling the acoustic signals received by the all-around ultrasonic transducer;
a data processing and imaging module coupled to the signal acquisition module, the data processing and imaging module including a memory and a processor, wherein the processor implements a method when executing a computer program stored in the memory, the method comprising:
scanning human organs based on a fully-enclosed ultrasonic transducer, and measuring the transit time of sound waves;
step two, reconstructing initial sound velocity distribution of the target;
step three, B ultrasonic imaging is carried out based on the sound velocity distribution of the target;
step four, extracting organ tissue boundary information from the B ultrasonic image;
step five, iteratively updating the propagation path of the sound wave;
step six, iteratively updating the sound velocity distribution of the target;
wherein:
the first step comprises the following steps: causing a fully-enclosed ultrasound transducer to transmit and receive ultrasound waves and extracting a transit time from the recorded ultrasound signals based on a "threshold-phase" method;
the second step comprises the following steps: the propagation path of the recording sound in the organ tissue is a straight line from the transmitting array element to the receiving array element, and is L(0)Represents;
recording the slowness of sound as s; and is
Solving formula L based on SART algorithm(0)s(0)Get an initial solution s of s ═ t(0)；
The third step comprises: based on the full-enclosure ultrasonic transducer, digital beam synthesis is carried out in a mode of electronic focusing and delay summation, and B ultrasonic imaging is carried out on human organs;
in the process of delay summation, calculating the emission delay time of each emission array element based on the sound slowness distribution of the organ tissues reconstructed in the second step or the fifth step;
the fourth step comprises the following steps: performing convolution and operation on each pixel point in the B ultrasonic imaging result by using a Sobel operator, and then selecting a threshold value to extract an edge; after the filtering operation is carried out, generating the boundary of the organ tissues;
the fifth step comprises the following steps: modeling refraction and reflection phenomena occurring on organ tissue boundaries during acoustic propagation based on acoustic slowness s(r)Calculating the sound arrival time of the ultrasonic wave reaching each point of the organ tissue boundary, and iteratively updating the propagation path of the ultrasonic wave in the organ tissue to obtain L(r+1)；
The sixth step comprises: iteratively updating sound velocity distribution of target, and solving L based on SART algorithm(r+1)s(r+1)Obtaining the updated sound slowness s as t(r+1)。
2. The system of claim 1, wherein the fully-enclosed ultrasound transducer comprises a ring-shaped ultrasound probe comprising a plurality of array elements, the system further comprising a lifting device, wherein the lifting device comprises:
the supporting frame is used for supporting the whole fully-enclosed ultrasonic transducer and the attached scanning equipment;
an inner cavity filled with an ultrasound couplant to receive a body organ; and
and the motor is used for driving the all-around ultrasonic transducer to perform linear motion along the normal direction of the reference plane of the human organ.
3. The system of claim 2, further comprising a control module for controlling the emission and reception of sound waves of the fully-enclosed ultrasonic transducer by switching an electronic switch on one hand and adjusting the height of the transducer by controlling the lifting device on the other hand to realize three-dimensional scanning of the human organ.
4. A computer-readable storage medium having a computer program stored thereon, the computer program when executed by a processor implementing a method comprising:
scanning human organs based on a fully-enclosed ultrasonic transducer, and measuring the transit time of sound waves;
step two, reconstructing initial sound velocity distribution of the target;
step three, B ultrasonic imaging is carried out based on the sound velocity distribution of the target;
step four, extracting organ tissue boundary information from the B ultrasonic image;
step five, iteratively updating the propagation path of the sound wave;
step six, iteratively updating the sound velocity distribution of the target;
wherein:
the first step comprises the following steps: causing a fully-enclosed ultrasound transducer to transmit and receive ultrasound waves and extracting a transit time from the recorded ultrasound signals based on a "threshold-phase" method;
the second step comprises the following steps: the propagation path of the recording sound in the organ tissue is a straight line from the transmitting array element to the receiving array element, and is L(0)Represents;
recording the slowness of sound as s; and is
Solving formula L based on SART algorithm(0)s(0)Get an initial solution s of s ═ t(0)；
The third step comprises: based on the full-enclosure ultrasonic transducer, digital beam synthesis is carried out in a mode of electronic focusing and delay summation, and B ultrasonic imaging is carried out on human organs;
in the process of delay summation, calculating the emission delay time of each emission array element based on the sound slowness distribution of the organ tissues reconstructed in the second step or the fifth step;
the fourth step comprises the following steps: performing convolution and operation on each pixel point in the B ultrasonic imaging result by using a Sobel operator, and then selecting a threshold value to extract an edge; after the filtering operation is carried out, generating the boundary of the organ tissues;
the fifth step comprises the following steps: modeling refraction and reflection phenomena occurring on organ tissue boundaries during acoustic propagation based on acoustic slowness s(r)Calculating the sound arrival time of the ultrasonic wave reaching each point of the organ tissue boundary, and iteratingPropagation path of new ultrasonic wave in organ tissue, obtaining L(r+1)；
The sixth step comprises: iteratively updating sound velocity distribution of target, and solving L based on SART algorithm(r+1)s(r+1)Obtaining the updated sound slowness s as t(r+1)。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201911012719.6A CN110772281B (en) | 2019-10-23 | 2019-10-23 | Ultrasonic CT imaging system based on improved ray tracing method |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201911012719.6A CN110772281B (en) | 2019-10-23 | 2019-10-23 | Ultrasonic CT imaging system based on improved ray tracing method |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN110772281A CN110772281A (en) | 2020-02-11 |
| CN110772281B true CN110772281B (en) | 2022-03-22 |
Family
ID=69386703
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201911012719.6A Active CN110772281B (en) | 2019-10-23 | 2019-10-23 | Ultrasonic CT imaging system based on improved ray tracing method |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN110772281B (en) |
Families Citing this family (9)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111214213B (en) * | 2020-02-13 | 2022-11-11 | 南京科技职业学院 | Photoacoustic tomography method suitable for medium with nonuniform sound velocity |
| CN112754527B (en) * | 2020-12-28 | 2023-10-20 | 沈阳工业大学 | Data processing method for low-frequency ultrasonic thoracic imaging |
| CN113109446B (en) * | 2021-04-15 | 2022-11-29 | 复旦大学 | Ultrasonic tomography method |
| GB2612093B (en) * | 2021-10-21 | 2023-11-08 | Darkvision Tech Inc | Ultrasonic inspection of complex surfaces |
| CN114847983B (en) * | 2022-04-14 | 2024-05-14 | 华中科技大学 | Breast ultrasound CT (computed tomography) reflection imaging parameter self-adaptive compensation method based on pre-scanning |
| CN114886469B (en) * | 2022-05-11 | 2024-06-14 | 中国科学院声学研究所 | Array element positioning method and device of ultrasonic CT array probe |
| CN115100164B (en) * | 2022-06-29 | 2024-09-13 | 华中科技大学 | Ultrasonic CT transition time automatic extraction method based on edge detection |
| CN115616084A (en) * | 2022-11-10 | 2023-01-17 | 浙江衡玖医疗器械有限责任公司 | Rapid simulation method for large-scale three-dimensional ultrasonic array data |
| CN117158911B (en) * | 2023-10-25 | 2024-01-23 | 杭州励影光电成像有限责任公司 | Multi-sound-velocity self-adaptive photoacoustic tomography image reconstruction method |
Citations (9)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101917908A (en) * | 2007-12-28 | 2010-12-15 | 博莱科瑞士股份有限公司 | Quantification analisys of immobilized contrast agent in medical imaging applications |
| CN104408398A (en) * | 2014-10-21 | 2015-03-11 | 无锡海斯凯尔医学技术有限公司 | Liver boundary identification method and system |
| CN105997153A (en) * | 2016-07-15 | 2016-10-12 | 华中科技大学 | Imaging method of ultrasonic CT |
| CN108171696A (en) * | 2017-12-29 | 2018-06-15 | 深圳开立生物医疗科技股份有限公司 | A kind of placenta detection method and device |
| JP2018134199A (en) * | 2017-02-21 | 2018-08-30 | 株式会社日立製作所 | Medical imaging apparatus, image processing method, and program |
| CN108490079A (en) * | 2018-03-19 | 2018-09-04 | 哈尔滨工业大学 | A kind of beam-forming method based on ultrasonic transducer |
| CN109646089A (en) * | 2019-01-15 | 2019-04-19 | 浙江大学 | A kind of spine and spinal cord body puncture based on multi-mode medical blending image enters waypoint intelligent positioning system and method |
| CN110051387A (en) * | 2019-04-11 | 2019-07-26 | 华中科技大学 | A kind of ultrasound computed tomography image rebuilding method and system based on ray theory |
| CN110353729A (en) * | 2019-07-30 | 2019-10-22 | 北京航空航天大学 | A kind of sound wave transition time detection method based on two-way shot and long term memory network |
Family Cites Families (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US9375194B2 (en) * | 2012-05-28 | 2016-06-28 | Doron Kwiat | Real-time localization of an interventional tool |
| US10231704B2 (en) * | 2013-12-20 | 2019-03-19 | Raghu Raghavan | Method for acquiring ultrasonic data |
| US20160135775A1 (en) * | 2014-11-17 | 2016-05-19 | Wisconsin Alumni Research Foundation | System And Method For Time-Resolved, Three-Dimensional Angiography With Physiological Information |
-
2019
- 2019-10-23 CN CN201911012719.6A patent/CN110772281B/en active Active
Patent Citations (9)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101917908A (en) * | 2007-12-28 | 2010-12-15 | 博莱科瑞士股份有限公司 | Quantification analisys of immobilized contrast agent in medical imaging applications |
| CN104408398A (en) * | 2014-10-21 | 2015-03-11 | 无锡海斯凯尔医学技术有限公司 | Liver boundary identification method and system |
| CN105997153A (en) * | 2016-07-15 | 2016-10-12 | 华中科技大学 | Imaging method of ultrasonic CT |
| JP2018134199A (en) * | 2017-02-21 | 2018-08-30 | 株式会社日立製作所 | Medical imaging apparatus, image processing method, and program |
| CN108171696A (en) * | 2017-12-29 | 2018-06-15 | 深圳开立生物医疗科技股份有限公司 | A kind of placenta detection method and device |
| CN108490079A (en) * | 2018-03-19 | 2018-09-04 | 哈尔滨工业大学 | A kind of beam-forming method based on ultrasonic transducer |
| CN109646089A (en) * | 2019-01-15 | 2019-04-19 | 浙江大学 | A kind of spine and spinal cord body puncture based on multi-mode medical blending image enters waypoint intelligent positioning system and method |
| CN110051387A (en) * | 2019-04-11 | 2019-07-26 | 华中科技大学 | A kind of ultrasound computed tomography image rebuilding method and system based on ray theory |
| CN110353729A (en) * | 2019-07-30 | 2019-10-22 | 北京航空航天大学 | A kind of sound wave transition time detection method based on two-way shot and long term memory network |
Non-Patent Citations (1)
| Title |
|---|
| 利用秩和条件数分析透射式超声CT的数据完备性;沈毅;《声学学报》;20120331;第37卷(第2期);第181-187页 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN110772281A (en) | 2020-02-11 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN110772281B (en) | Ultrasonic CT imaging system based on improved ray tracing method | |
| JP7204493B2 (en) | Mapping within a body cavity using an ultrasound array distributed over a basket catheter | |
| CN112771374B (en) | Training-based nonlinear mapping image reconstruction method | |
| CN107613881B (en) | Method and system for correcting fat-induced aberrations | |
| CN110785126B (en) | Method for ultrasound system independent attenuation coefficient estimation | |
| CN111819467B (en) | Method and apparatus for estimating wave propagation and scattering parameters | |
| US10497284B2 (en) | Systems and methods of ultrasound simulation | |
| JP2021501656A (en) | Intelligent ultrasound system to detect image artifacts | |
| WO2013116809A1 (en) | Ultrasound waveform tomography with tv regularization | |
| CN109589131B (en) | Ultrasonic method and ultrasonic system for automatically setting Doppler imaging mode parameters in real time | |
| US20130041261A1 (en) | Method and system for multi-grid tomographic inversion tissue imaging | |
| CN114972567B (en) | Medical ultrasonic CT multi-parameter image reconstruction method based on wave equation | |
| JP2021530303A (en) | Ultrasound imaging with deep learning and related devices, systems, and methods | |
| CA3046839A1 (en) | Method of, and apparatus for, non-invasive medical imaging using waveform inversion | |
| Perez-Liva et al. | Speed of sound ultrasound transmission tomography image reconstruction based on Bézier curves | |
| WO2017163389A1 (en) | Ultrasonic ct apparatus and ultrasonic imaging method | |
| CN117158911B (en) | Multi-sound-velocity self-adaptive photoacoustic tomography image reconstruction method | |
| KR20140137037A (en) | ultrasonic image processing apparatus and method | |
| KR102490019B1 (en) | Method and apparatus for quantitative imaging using ultrasound data | |
| JP7356504B2 (en) | Ultrasonic estimation of nonlinear bulk elasticity of materials | |
| JP2020502517A (en) | Method for acquiring a signal by ultrasonic search, corresponding computer program and ultrasonic search device | |
| CN111223157A (en) | Ultrasonic CT sound velocity imaging method based on depth residual error network | |
| JP5823184B2 (en) | Ultrasonic diagnostic apparatus, medical image processing apparatus, and medical image processing program | |
| KR102704209B1 (en) | Ultrasound image apparatus and method for controlling thereof | |
| CN114098795B (en) | System and method for generating ultrasound probe guidance instructions |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |