An object boundary extraction method applicable to ultrasonic and photoacoustic tomography

Abstract

The present invention provides a method for extracting the boundary of an object applicable to ultrasonic and photoacoustic tomography. The present invention proposes a method for directly extracting the boundary of an imaging object from ultrasonic / photoacoustic channel signals, deducing the possible spatial positions of the object boundary based on the transit time of the first-arrival echo, and accurately restoring the object boundary position by using methods such as spatial outlier rejection, smoothing, and interpolation, which is applicable to ultrasonic and photoacoustic tomography. Compared with the method of digital image processing, the method provided by the present invention is not affected by the quality of the reconstructed image and can quickly provide accurate boundary position information without reconstructing the image.

Description

Technical Field

[0001] The present invention belongs to the technical fields of ultrasonic tomography and photoacoustic tomography. Background Art

[0002] Ultrasonic tomography is a new type of ultrasonic imaging method. It not only retains the advantages of traditional ultrasonic imaging such as non-invasive, non-radiative, and good safety, but also, due to its unique tomographic scanning imaging method, can provide echo images with high spatial resolution and quantitative images reflecting the sound velocity and attenuation of the imaged object. A typical ultrasonic tomography system uses an annular ultrasonic transducer array to detect ultrasonic signals, emits divergent waves in the way of synthetic aperture to image an object, the elements of each channel emit ultrasonic waves in turn, and all channels receive ultrasonic signals together. This imaging method has a large number of acquisitions, a large amount of data, and a slow imaging speed. The method of plane wave emission can scan the entire imaging area in a relatively short time, greatly reducing the acquisition time, but the imaging quality is poor. Plane wave imaging can be achieved by using a polygonal linear array probe instead of an annular array, or by applying a time delay to the annular array to make it emit plane waves. Photoacoustic tomography combines the advantages of optical and acoustic imaging. By using the light absorption differences of biological tissues, it can provide structural and functional information inside the tissues, and also often uses an annular array to detect acoustic signals.

[0003] In ultrasonic and photoacoustic tomography technologies, obtaining the boundary information of the imaged object a priori helps to reduce the number of unknown variables in the image reconstruction area, improve the quality and speed of image reconstruction, and can also be used for foreground and background segmentation of post-processed images. At present, for boundary extraction in ultrasonic and photoacoustic tomography, digital image processing methods (such as algorithms like Snakes, active contours, and Grabcut) are mainly used to extract the boundaries of the reconstructed images, and accurate boundary information cannot be obtained before image reconstruction, thus unable to improve the image reconstruction process. In addition, due to differences in imaging devices and reconstruction parameters, the results of image boundary segmentation highly depend on the quality of the reconstructed images and are easily affected by background noise and artifacts. Summary of the Invention

[0004] The technical problem to be solved by the present invention is that the boundary extraction methods in ultrasonic and photoacoustic tomography cannot obtain accurate boundary information before image reconstruction, and the results of image boundary segmentation highly depend on the quality of the reconstructed images and are easily affected by background noise and artifacts.

[0005] To solve the above technical problem, the technical solution of the present invention is to provide an object boundary extraction method applicable to ultrasonic and photoacoustic tomography, which is characterized by including the following steps:

[0006] Step 1: Image the imaging target according to the requirements of ultrasonic tomography or photoacoustic tomography, and receive ultrasonic signals or photoacoustic signals through an ultrasonic transducer array; select the acquisition signals of specific transducer receiving channels for subsequent boundary extraction according to the imaging sequences generated by different imaging methods, where:

[0007] For ultrasonic tomography based on divergent wave emission: select appropriate emission-reception pair signals for subsequent boundary extraction;

[0008] For photoacoustic tomography: select the acquisition signals of all transducer receiving channels for subsequent boundary extraction;

[0009] For ultrasonic tomography based on plane wave emission: emit a plane wave at a certain angle, and select the acquisition signals of the transducer emission channel and the transducer receiving channels near it for subsequent boundary extraction;

[0010] Step 2: Extract the transit time of the first-arrival echo from the acquisition signals of the specific transducer receiving channels selected in Step 1;

[0011] Step 3: Let r rx =(x rx , z rx ) be the spatial coordinates of a certain transducer receiving channel, r=(x, z) be the reflection point on the imaging target, and c be the sound speed of the coupling medium. Then, the transit time of the first-arrival echo generated by the boundary is expressed by the following formula (1), (2), or (3):

[0012] For ultrasonic tomography based on divergent wave emission:

[0013]

[0014] In formula (1), r tx is the spatial coordinate of the divergent wave source, and τ SA (r tx , r rx , r) is the transit time of the first-arrival echo generated by the boundary;

[0015] For photoacoustic tomography:

[0016]

[0017] In formula (2), τ PA (r rx , r) is the transit time of the first-arrival echo generated by the boundary;

[0018] For ultrasonic tomography based on plane wave emission:

[0019]

[0020] In Equation (3), α is the plane wave emission angle, that is, the angle between the plane wave and the surface of the linear array ultrasonic probe. Here, the linear array ultrasonic probe is a real linear array ultrasonic probe or a virtual linear array ultrasonic probe obtained by superimposing and delaying ultrasonic arrays of other shapes; L is the aperture size of the linear array ultrasonic probe; τ PW (α, r rx , r) is the transit time of the first arrival echo generated by the boundary;

[0021] Among them, arranging Equation (1), (2) or (3) into the form of an ellipse equation represents all potential positions of the boundary points deduced from the transit time of the first arrival echo, as shown in the following Equation (4):

[0022] Ax 2 +Bxz + Cz 2 +Dx + Ez + F = 0 (4)

[0023] In Equation (4), A - F are constants;

[0024] For the ellipse equations corresponding to two adjacent transducer receiving channels, calculate the common tangent and the corresponding tangent points, and take the tangent points as the approximation of the boundary points, which are defined as the initial boundary points;

[0025] Traverse the specific transducer receiving channels selected in Step 1 in sequence, calculate all the initial boundary points, and obtain the initial boundary point set;

[0026] Step 4: Post-process the initial boundary point set to obtain the final boundary detection result, including the following steps:

[0027] Step 401: Assign the initial boundary points in the initial boundary point set to M sector regions evenly distributed at a certain angle in the polar coordinate system;

[0028] Step 402: In each sector region, perform spatial outlier rejection on the initial boundary points, and only retain the initial boundary point with the shortest distance from the origin of the polar coordinates;

[0029] Step 403: Perform smoothing and denoising processing on the discrete initial boundary points after outlier rejection;

[0030] Step 404: Interpolate the smoothed discrete initial boundary points, and all the boundary points obtained are the final boundary detection result.

[0031] Preferably, in Step 2, for any acquisition signal with a sampling length of N, the following steps are used to extract the transit time TOF of the first arrival echo of the acquisition signal:

[0032] Step 201: Calculate the AIC value of each sampling point. The calculation formula for the AIC value AIC(k) of the sampling point at the k-th moment is shown in the following formula (5):

[0033]

[0034] In formula (8), and are the variances of the two segments of the acquired signals within the time periods [1, k] and [k + 1, N] respectively;

[0035] Step 202: Among all the sampling points, calculate the time of flight TOF of the first arriving echo based on the moment corresponding to the sampling point with the minimum AIC value, as shown in the following formula (6):

[0036] TOF = T s *argmin k {AIC(k)} + t0 (6)

[0037] In formula (9), T s is the sampling interval, and t0 is the starting time of the acquired signal with the current sampling length of N. Before detection, the acquired signals of the channels with relatively low signal-to-noise ratio can be discarded, and the acquired signals of other channels with high signal-to-noise ratio are used to interpolate the missing channels' acquired signals. To further reduce the interference of noise, the acquired signals can also be smoothed to improve the accuracy of the subsequent steps.

[0038] Preferably, in step 3, the ultrasonic tomography imaging system based on divergent wave emission is composed of the following four methods:

[0039] Method 1): Divergent waves are directly generated by a single or multiple array elements of a circular array ultrasonic transducer. In the present invention, the circular array includes a 360-degree closed circular array and also includes an open circular array such as 90 degrees, and the same applies hereinafter and will not be elaborated further;

[0040] Method 2): Delays are superimposed on some array elements of the circular ultrasonic array to obtain divergent waves generated by a virtual source;

[0041] Method 3): Divergent waves are directly generated by a single array element of a linear array ultrasonic transducer, or delays are superimposed on multiple array elements of the linear array ultrasonic transducer to generate divergent waves, and the probe or the imaging target is rotated for tomographic scanning; In the present invention, the linear array ultrasonic transducer uses a flat probe, and the same applies hereinafter and will not be elaborated further;

[0042] Method 4): Divergent waves are directly generated by a single or multiple array elements of a linear convex array ultrasonic transducer or through superimposing delays, and the probe or the imaging target is rotated for tomographic scanning.

[0043] Preferably, in step 3, the ultrasonic tomography system based on plane wave emission is composed of the following three methods:

[0044] Method 1) A polygon array is formed by multiple conventional linear array ultrasonic probes;

[0045] Method 2) Tomographic scanning is performed by rotating a single or multiple conventional linear array ultrasonic probes or rotating the imaging target;

[0046] Method 3) Delay is superimposed on the annular ultrasonic array to make it emit a plane wave, then:

[0047] The formula (3) is calculated in the coordinate system of the currently emitted linear array ultrasonic probe or the virtual linear array ultrasonic probe obtained by superimposing the delay of the annular ultrasonic array.

[0048] Preferably, in step 3, for ultrasonic tomography based on plane wave emission, the tangent point is calculated in the coordinate system of the emitted linear array probe, and then the coordinates are transformed to a unified global coordinate system.

[0049] Preferably, in step 3, for ultrasonic tomography based on divergent wave emission, A-F in formula (4) is obtained from r tx 、r rx 、c and τ SA Calculated.

[0050] Preferably, in step 3, for ultrasonic tomography based on divergent wave emission, when the transducer for receiving signals is also the transmitting channel, then: r tx =r rx ,A=C=1,B=0,D=-2x rx ,E=-2z rx , Formula (4) is simplified to:

[0051]

[0052] In formula (7), τ SA =τ SA (r tx =r rx ,r rx ,r).

[0053] Preferably, in step 3, for photoacoustic tomography, A=C=1,B=0,D=-2x rx ,E=-2z rx , Formula (4) is expressed as:

[0054]

[0055] In formula (8), τ PA= τ PA (r rx , r).

[0056] Preferably, in step 3, for ultrasonic tomography based on plane wave emission, A - F in formula (4) are calculated from α, L, r rx , c and τ PW .

[0057] Preferably, in step 3, for ultrasonic tomography based on plane wave emission, when the plane wave angle α = 0°, A = 1, B = C = 0, D = -2x rx , E = 2(cτ PW -z rx ), Formula (4) is simplified to;

[0058]

[0059] In formula (9), τ PW = τ PW (0°, r rx , r).

[0060] Preferably, in step 4, when performing spatial outlier rejection on the initial boundary points, the fan-shaped region is deflected by the angle θ to obtain other fan-shaped regions; subsequently, the initial boundary points with the shortest distance to the origin in all the deflected fan-shaped regions are retained to achieve outlier rejection.

[0061] Preferably, in step 4, during the post-processing, according to the distribution of the outliers, the parameters M and θ are adjusted to retain the detection results closer to the boundary of the imaging target.

[0062] The present invention proposes a method for directly extracting the boundary of an imaging object from ultrasonic / photoacoustic channel signals, deducing the possible spatial positions of the object boundary based on the transit time of the first-arrival echo, and accurately restoring the object boundary position using methods such as spatial outlier rejection, smoothing, and interpolation, which is applicable to ultrasonic and photoacoustic tomography. Compared with the method of digital image processing, the method provided by the present invention is not affected by the quality of the reconstructed image and can quickly provide accurate boundary position information without reconstructing the image. BRIEF DESCRIPTION OF THE DRAWINGS

[0063] Figure 1 is the flowchart of the boundary extraction algorithm;

[0064] Figure 2 is the schematic diagram of the coordinate system for ultrasonic tomography based on plane wave emission;

[0065] Figure 3 is the schematic diagram of outlier rejection of boundary points;

[0066] Figure 4 Schematically shows the transit time of the first-arrival echo of finger ultrasound tomography;

[0067] Figure 5 Schematically shows a finger boundary detector for ultrasound tomography based on divergent wave emission, where: (a) the initial boundary points after calculating the tangent points for two adjacent receiving channels; (b) the boundary points after spatial outlier rejection; (c) the spatial positions of the discrete boundary points after smoothing out the outliers; (d) interpolating the spatial positions of the discrete boundary points to obtain the final boundary detection result;

[0068] Figure 6 Schematically shows the transit time of the first-arrival echo of finger photoacoustic tomography;

[0069] Figure 7 Schematically shows finger boundary detection for photoacoustic tomography, where: (a) the initial boundary points after calculating the tangent points for two adjacent receiving channels; (b) the final boundary detection result after post-processing;

[0070] Figure 8 Schematically shows breast boundary detection for plane-wave ultrasound tomography, where: (a) eight linear arrays form an octagonal ultrasound tomography system; (b) the sound speed distribution map of the breast model set in the simulation; (c) the initial boundary points after calculating the tangent points for two adjacent receiving channels; (d) the final boundary detection result after smoothing and interpolation processing. Specific embodiments

[0071] The present invention will be further described below in conjunction with specific embodiments. It should be understood that these embodiments are only used to illustrate the present invention and not to limit the scope of the present invention. In addition, it should be understood that after reading the content taught by the present invention, those skilled in the art can make various changes or modifications to the present invention, and these equivalent forms also fall within the scope defined by the appended claims of this application.

[0072] As Figure 1 shown, an object boundary extraction method applicable to ultrasound and photoacoustic tomography provided by the present invention includes the following steps:

[0073] Step 1: Data acquisition

[0074] Image the imaging target according to the imaging requirements of ultrasound or photoacoustics, and receive ultrasound or photoacoustic signals through an ultrasound transducer array. According to different imaging sequences, select the acquisition signals of specific transducer receiving channels for subsequent boundary extraction, where:

[0075] For the ultrasound tomography sequence based on divergent wave emission: Select appropriate transmit-receive pair signals for subsequent boundary extraction, taking the acquisition signals of the transducer receiving channels that are both transmitters and receivers as an example;

[0076] For the photoacoustic tomography sequence: The acquisition signals of all receiving channels are selected for subsequent boundary extraction;

[0077] For the ultrasound tomography sequence based on plane wave emission: A plane wave at a certain angle (taking the 0 ° plane wave as an example) is emitted, and the acquisition signals of the transducer emission channel and the transducer receiving channels near it are selected for subsequent boundary extraction.

[0078] Step 2: Extraction of the time-of-flight of the first-arrival echo

[0079] Methods such as the Akaike information criterion (AIC) are used to extract the time-of-flight (TOF) of the first-arrival echo from the acquisition signals of the specific receiving channels selected in Step 1. For any acquisition signal with a sampling length of N, there is:

[0080] Calculate the AIC value of each sampling point. Among them, the calculation formula for the AIC value AIC(k) of the sampling point at the k-th moment is shown in the following formula (1):

[0081]

[0082] In formula (1), and are the variances of two segments of acquisition signals within the time periods [1, k] and [k + 1, N] respectively;

[0083] Among all sampling points, the time corresponding to the sampling point with the minimum AIC value is calculated as the time-of-flight TOF of the first-arrival echo, as shown in the following formula (2):

[0084] TOF = T s *argmin k {AIC(k)} + t0 (2)

[0085] In formula (2), T s is the sampling interval, and t0 is the starting time of the acquisition signal with the current sampling length of N.

[0086] The time-of-flight of the first-arrival echo can also be calculated using methods such as first crossing zero extraction, cross-correlation, and first echo peak extraction, which are mainly divided into the following four categories:

[0087] 1. The method based on template matching, such as cross-correlation; 2. The method based on evaluation functions: Using evaluation functions such as the Akaike information criterion, energy ratio, and entropy as the judgment criteria for different features near the time-of-flight point; 3. The method based on neural networks. 4. The method based on the characteristic points of the signal to be measured: such as first crossing zero extraction, first echo peak extraction, first significant peak extraction, etc.

[0088] Due to the presence of noise and crosstalk in the original signal, the detection result of the transit time is easily affected. Before detection, the acquisition signals of channels with low signal-to-noise ratio can be discarded, and the acquisition signals of other channels with high signal-to-noise ratio are used to interpolate the acquisition signals of the missing channels. To further reduce the interference of noise, the acquisition signals can also be smoothed to improve the accuracy of subsequent steps.

[0089] Step 3: Initial boundary point extraction

[0090] Establish the connection between the boundary spatial position and the transit time TOF of the first-arrival echo. Assume r rx =(x rx , z rx ) is the spatial coordinate of a certain receiving channel, r = (x, z) is the reflection point on the imaging object, and c is the sound speed of water. Then, the transit time of the first-arrival echo generated by the boundary can be expressed by the following formula (3), (4), or (5):

[0091] For ultrasonic tomography based on divergent wave emission:

[0092]

[0093] In formula (3), r tx is the spatial coordinate of the receiving channel, τ SA (r tx = r rx , r rx , r) is the transit time of the first-arrival echo generated by the boundary. Among them, the ultrasonic tomography system based on divergent wave emission is composed of the following four methods:

[0094] Method 1) Divergent waves are directly generated by a single or multiple elements of a circular array ultrasonic transducer. It should be noted that;

[0095] Method 2) Delays are superimposed on some elements of the circular ultrasonic array to obtain divergent waves generated by a virtual source;

[0096] Method 3) Divergent waves are directly generated by a single element of a linear array ultrasonic transducer, or delays are superimposed on multiple elements of the linear array ultrasonic transducer to generate divergent waves, and the probe or the imaging target is rotated for tomographic scanning;

[0097] Method 4) Divergent waves are directly generated by a single or multiple elements of a linear convex array ultrasonic transducer or by superimposing delays, and the probe or the imaging target is rotated for tomographic scanning. Formula (3) holds for all the above four methods.

[0098] For photoacoustic tomography:

[0099]

[0100] In Equation (4), τ PA (r rx , r) is the transit time of the first-arrival echo generated by the boundary;

[0101] For ultrasonic tomography based on plane-wave emission:

[0102]

[0103] In Equation (5), τ PW (0°, r rx , r) is the transit time of the first-arrival echo generated by the boundary, where the ultrasonic tomography system based on plane-wave emission can be composed of the following three methods:

[0104] Method 1) A polygon array is composed of multiple conventional linear array ultrasonic probes;

[0105] Method 2) Tomographic scanning is performed by rotating a single or multiple conventional linear array ultrasonic probes or rotating the imaging target;

[0106] Method 3) Delays are added to the annular ultrasonic array to make it emit plane waves.

[0107] Equation (5) holds for all of the above three methods, and it is required to calculate in the coordinate system (local coordinate system) where the currently emitted linear array ultrasonic probe (or the virtual linear array ultrasonic probe obtained by adding delays to the annular ultrasonic array) is located, as Figure 2 shown.

[0108] Equations (3) - (5) are arranged into the form of an ellipse equation, representing all potential positions of the boundary points deduced from the transit time of the first-arrival echo, where:

[0109] For ultrasonic tomography based on divergent-wave emission:

[0110]

[0111] In Equation (6), τ SA = τ SA (r tx = r rx , r rx , r);

[0112] For photoacoustic tomography:

[0113]

[0114] In Equation (7), τ PA = τ PA (r rx , r);

[0115] For ultrasonic tomography based on plane wave emission:

[0116]

[0117] In Equation (8), τ PW = τ PW (0°, r rx , r).

[0118] For the elliptic equations shown in Equations (6), (7), or (8) corresponding to two adjacent receiving channels, calculate the common tangent and the corresponding tangent points, and use the tangent points as the approximation of the boundary points.

[0119] Traverse the specific receiving channels selected in Step 1 in sequence, calculate the tangent points. All the calculated tangent points form the initial boundary point set, and each point in the initial boundary point set is further defined as an initial boundary point. For ultrasonic tomography based on plane wave emission, the tangent point calculation needs to be performed in the coordinate system (local coordinate system) where the linear array probe for emission is located. After the tangent point calculation is completed, coordinate transformation is performed to the unified global coordinate system.

[0120] Step Four: Post-processing

[0121] Allocate the detected initial boundary points into M sector regions evenly distributed at a certain angle in the polar coordinate system. In each sector region, perform spatial outlier rejection on the initial boundary points, and only retain the initial boundary point with the shortest distance from the origin of the polar coordinates. The number of retained boundary points is determined by the set number M of sector regions. Among them, in order to improve the effect of outlier rejection, the angle of each sector region can be increased by θ without changing the number M of sector regions. For example, Figure 3 in it, deflect the angle of the black sector region by θ to obtain two other sector regions. Retain the initial boundary points with the shortest distance from the origin in all the sector regions after the angle deflection to achieve outlier rejection.

[0122] Perform smoothing and denoising processing on the discrete initial boundary points after outlier rejection.

[0123] Finally, interpolate the smoothed discrete initial boundary points, and all the obtained boundary points are the final boundary detection results.

[0124] During the post-processing, according to the distribution of outliers, the parameters M and θ can be adjusted to retain the detection results closer to the object boundary.

[0125] Example 1: Ultrasonic tomography based on divergent wave emission

[0126] An object boundary extraction method applicable to ultrasonic and photoacoustic tomography disclosed in this example includes the following steps:

[0127] Step 1: Data acquisition

[0128] Use a 512-channel annular ultrasound array with a radius of 40 mm to perform ultrasonic imaging on the finger. The central frequency is 5 MHz and the sampling frequency is 20 MHz. Using the synthetic aperture method, each channel element sequentially emits divergent waves and all channels receive. Use the data of 128 channels centered on the transmitting element to reconstruct the ultrasonic echo image, and only use the data received by the transmitting element to reconstruct the finger boundary.

[0129] Step 2: Extraction of the transit time of the first-arrival echo

[0130] Use the AIC method to extract the transit time of the first-arrival echo, and the results are as Figure 4 shown.

[0131] Step 3: Extraction of initial boundary points

[0132] According to formula (6), calculate the common tangents and tangent points pairwise for the elliptic equations between adjacent two receiving channels in sequence to obtain the initial boundary point set, as shown in Figure 5 (a).

[0133] Step 4: Post-processing

[0134] Perform spatial outlier rejection on the initial boundary point set, where the number M of fan-shaped regions is set to 128 and the angular offset θ of the fan-shaped regions is set to 1°. The results are as shown in Figure 5 (b), where all outliers are rejected and the remaining boundary points are closer to the original boundary. To reduce the influence of noise, perform spatial coordinate smoothing on the discrete boundary points after outlier rejection, and finally obtain the final boundary detection result through spatial coordinate interpolation, and the results are as shown in Figure 5 (c) and (d).

[0135] Example 2: Photoacoustic tomography

[0136] A method for extracting the boundary of an object applicable to ultrasonic and photoacoustic tomography disclosed in this example includes the following steps:

[0137] Step 1: Data acquisition

[0138] Use a 128-channel annular ultrasound array with a radius of 30 mm to perform photoacoustic imaging on the finger. The central frequency is 2.5 MHz and the sampling rate is 40 MHz. The laser emission pulse repetition frequency is 10 Hz and the wavelength is 720 nm. Interpolate the 128-channel data to 512 channels to improve the quality of photoacoustic image reconstruction and the distribution density of initial boundary points.

[0139] Step 2: Extraction of the transit time of the first-arrival wave

[0140] The transit time of the first-arrival echo is extracted using the AIC method, and the transit time is smoothed using a 7-point mean smoothing filter. The result is as Figure 6 shown.

[0141] Step 3: Initial boundary point extraction

[0142] According to formula (7), the common tangents and tangent points are calculated pairwise for the elliptic equations between adjacent two receiving channels in sequence, and an initial boundary point set is obtained, as Figure 7 (a) shown.

[0143] Step 4: Post-processing

[0144] Spatial outlier rejection is performed on the initial boundary point set, and smoothing and interpolation processing are performed on the discrete coordinates after outlier rejection to obtain the final boundary detection result, as Figure 7 (b) shown.

[0145] Embodiment 3: Ultrasonic tomography based on plane wave emission

[0146] An object boundary extraction method applicable to ultrasonic and photoacoustic tomography disclosed in this embodiment includes the following steps:

[0147] Step 1: Data acquisition

[0148] In the simulation software k-Wave, a numerical model of a breast slice is subjected to plane wave ultrasonic tomography using a 1024-channel octagonal ultrasonic array with a radius of 40 mm, and the center frequency is set to 4.5 MHz. Eight linear arrays surround the breast model, as Figure 8 (a) and Figure 8 (b) shown. Each array emits 5 plane waves in sequence, and the tilt angles of the emitted plane waves vary from -5° to 5° at intervals of 2.5°. Three linear arrays centered on the emission array receive echo signals for image reconstruction, as Figure 2 shown in. The echo signal received by the emission array when the 0° plane wave is emitted is used for boundary detection.

[0149] Step 2: First-arrival wave transit time extraction

[0150] The transit time of the first-arrival echo is extracted using the AIC method, and the transit time is smoothed using a 5-point mean smoothing filter.

[0151] Step 3: Initial boundary point extraction

[0152] According to formula (8), the common tangents and tangent points are calculated pairwise for the elliptic equations between adjacent two receiving channels in sequence, and an initial boundary point set is obtained, as Figure 8 (c) shown.

[0153] Step 4: Post-processing

[0154] Spatial outlier rejection is performed on the initial boundary point set, and smoothing and interpolation are performed on the discrete coordinates after outlier rejection to obtain the final boundary detection result, as Figure 8 (d) shown.

Claims

1. An object boundary extraction method applicable to ultrasonic and photoacoustic tomography, characterized in that It includes the following steps: Step 1: Image the imaging target according to the requirements of ultrasonic tomography or photoacoustic tomography, and receive ultrasonic signals or photoacoustic signals through an ultrasonic transducer array; according to the imaging sequences generated by different imaging methods, select the acquisition signals of specific transducer receiving channels for subsequent boundary extraction, where: For ultrasonic tomography based on divergent wave emission: Select appropriate transmit-receive pair signals for subsequent boundary extraction; For photoacoustic tomography: Select the acquisition signals of all transducer receiving channels for subsequent boundary extraction; For ultrasonic tomography based on plane wave emission: For ultrasonic tomography based on plane wave emission: Emit a plane wave at a certain angle, and select the acquisition signals of the transducer transmitting channel and the transducer receiving channels near it for subsequent boundary extraction; Step 2: Extract the transit time of the first-arrival echo from the acquisition signals of the specific transducer receiving channels selected in Step 1; Step 3. Let r rx =(x rx , z rx ) be the spatial coordinates of a receiving channel of a transducer, r = (x, z) be the reflection point on the imaging target, and c be the sound speed of the coupling medium. Then, the transit time of the first-arrival echo generated by the boundary is expressed by the following formula (1), (2), or (3): For ultrasonic tomography based on divergent wave emission: In Equation (1), r tx is the spatial coordinate of the divergent wave source, and τ SA (r tx , r rx , r) is the transit time of the first arrival echo generated by the boundary; For photoacoustic tomography: In formula (2), τ PA (r rx , r) is the transit time of the first arriving echo generated by the boundary; For ultrasonic tomography based on plane wave emission: In Equation (3), α is the plane wave emission angle, that is, the angle between the plane wave and the surface of the linear array ultrasonic probe, where the linear array ultrasonic probe is a real linear array ultrasonic probe or a virtual linear array ultrasonic probe obtained by superimposing and delaying other shaped ultrasonic arrays; L is the aperture size of the linear array ultrasonic probe; τ PW (α, r rx , r) is the transit time of the first arrival echo generated by the boundary; Among them, arranging Equation (1), (2) or (3) into the form of an ellipse equation represents all potential positions of the boundary points deduced from the transit time of the first-arrival echo, as shown in the following Equation (4): Ax 2 + Bxz + Cz 2 + Dx + Ez + F = 0 (4) In Equation (4), A - F are constants; Traverse the specific transducer receiving channels selected in Step 1 in sequence, calculate all the initial boundary points, and obtain the initial boundary point set; Step 4: Post-process the initial boundary point set to obtain the final boundary detection result, including the following steps: Step 401: Assign the initial boundary points in the initial boundary point set to M fan-shaped regions evenly distributed at a certain angle in the polar coordinate system; Step 402: In each fan-shaped region, perform spatial outlier rejection on the initial boundary points, and only retain the initial boundary point with the shortest distance from the origin of the polar coordinates; Step 403: Perform smoothing and denoising processing on the discrete initial boundary points after outlier rejection; Step 404: Interpolate the smoothed discrete initial boundary points, and all the boundary points obtained are the final boundary detection results.

2. The object boundary extraction method applicable to ultrasonic and photoacoustic tomography according to claim 1, characterized in that In Step 2, for any acquisition signal with a sampling length of N, the following steps are used to extract the transit time TOF of the first-arrival echo of the acquisition signal: Step 201: Calculate the AIC value of each sampling point. Among them, the calculation formula for the AIC value AIC(k) of the sampling point at the k-th moment is shown in the following Equation (5): In formula (5), and are the variances of two segments of the acquired signals in time periods [1, k] and [k + 1, N], respectively; Step 202: Among all the sampling points, calculate the transit time TOF of the first-arrival echo based on the moment corresponding to the sampling point with the minimum AIC value, as shown in the following Equation (6): TOF = T s *argmin k {AIC(k)} + t0 (6) In Equation (6), T s is the sampling interval, and t0 is the starting time of the acquisition signal with the current sampling length of N; Before detection, discard the acquisition signals of channels with low signal-to-noise ratio, and use the acquisition signals of other channels with high signal-to-noise ratio to interpolate the acquisition signals of the missing channels; in order to further reduce the interference of noise, the acquisition signals are also smoothed to improve the accuracy of subsequent steps.

3. The object boundary extraction method applicable to ultrasonic and photoacoustic tomography according to claim 1, characterized in that, In Step 3, the ultrasonic tomography system based on divergent wave emission is composed of the following four methods: Method 1) Directly generate divergent waves by a single or multiple array elements of a circular array ultrasonic transducer; Method 2) Delay is superimposed on some elements of the annular ultrasonic array to obtain a divergent wave generated by a virtual source; Method 3) A divergent wave is directly generated by a single element of the linear array ultrasonic transducer, or a divergent wave is generated by superimposing time delays on multiple elements of the linear array ultrasonic transducer, and the probe or the imaging target is rotated for tomographic scanning; Method 4) A divergent wave is directly generated by a single or multiple elements of the linear convex array ultrasonic transducer or through superimposing time delays, and the probe or the imaging target is rotated for tomographic scanning.

4. A method for extracting the object boundary applicable to ultrasonic and photoacoustic tomography according to claim 1, characterized in that In step 3, the ultrasonic tomographic imaging system based on plane wave emission is composed of the following three methods: Method 1) A polygon array is formed by multiple conventional linear array ultrasonic probes; Method 2) Tomographic scanning is performed by rotating a single or multiple conventional linear array ultrasonic probes or rotating the imaging target; Method 3) Delay is superimposed on the annular ultrasonic array to make it emit a plane wave, then: Equation (3) is calculated in the coordinate system of the currently emitted linear array ultrasonic probe or the virtual linear array ultrasonic probe obtained by superimposing time delays on the annular ultrasonic array.

5. A method for extracting the object boundary applicable to ultrasonic and photoacoustic tomography according to claim 1, characterized in that, In step 3, for ultrasonic tomography based on plane wave emission, the tangent point is calculated in the coordinate system of the emitted linear array probe, and then the coordinates are transformed to a unified global coordinate system.

6. The object boundary extraction method applicable to ultrasonic and photoacoustic tomography according to claim 1, wherein In step 3, for ultrasonic tomography based on divergent wave emission, when the transducer for receiving signals is also the transmitting channel, there is: r tx = r rx , A = C = 1, B = 0, D = -2x rx , E = -2z rx , Formula (4) is simplified to: In formula (7), τ SA = τ SA (r tx = r rx , r rx , r).

7. The object boundary extraction method applicable to ultrasonic and photoacoustic tomography according to claim 1, characterized in that In step 3, for photoacoustic tomography, A = C = 1, B = 0, D = -2x rx , E = -2z rx , Equation (4) is expressed as: In Equation (8), τ PA = τ PA (r rx , r).

8. The object boundary extraction method applicable to ultrasonic and photoacoustic tomography according to claim 1, wherein For ultrasonic tomography based on plane wave emission, when the plane wave angle α = 0°, A = 1, B = C = 0, D = -2x rx , E = 2(cτ PW -z rx ), Formula (4) is simplified to; In Equation (9), τ PW = τ PW (0°, r rx , r).

9. The object boundary extraction method applicable to ultrasonic and photoacoustic tomography according to claim 1, wherein In step 4, when performing spatial outlier rejection on the initial boundary points, the fan-shaped region is deflected by an angle θ to obtain other fan-shaped regions; subsequently, the initial boundary points with the shortest distance from the origin in all the fan-shaped regions after the angle deflection are retained to achieve outlier rejection.

10. A method for extracting the object boundary applicable to ultrasonic and photoacoustic tomography according to claim 1, characterized in that, In step 4, during the post-processing, according to the distribution of the outliers, the parameters M and θ are adjusted to retain the detection results closer to the boundary of the imaging target.

Images