Submarine karst morphology detection method based on multi-frequency ultrasound

Measurement 155 (2020) 107532

Contents lists available at ScienceDirect

Measurement journal homepage: www.elsevier.com/locate/measurement

Submarine karst morphology detection method based on multi-frequency ultrasound Jin-chao Wang a, Chuan-ying Wang a, Zeng-qiang Han a, Yi-teng Wang a, Xiao-hua Huang b,⇑ a b

State Key Laboratory of Geomechanics and Geotechnical Engineering, Institute of Rock and Soil Mechanics, Chinese Academy of Sciences, Wuhan, Hubei 430071, China Guangxi Key Laboratory of Disaster Prevention and Engineering Safety, Guangxi University, Nanning, Guangxi 530004, China

a r t i c l e

i n f o

Article history: Received 26 May 2019 Received in revised form 18 December 2019 Accepted 20 January 2020 Available online 23 January 2020 Keywords: Submarine karst Multi-frequency ultrasound Working frequency Morphology detection Automatic matching

a b s t r a c t Karst is a geological phenomena in which water mainly dissolves soluble rocks, supplemented by mechanical actions, such as water erosion, potential erosion and collapse. During the design and construction of submarine tunnel, the fine detection of karst cave shape plays an important role, which is related to the construction process of karst area and the stability evaluation in the later stage. However, due to the irregularity and complexity of submarine karst morphology, it is difficult for the existing detection technology to realize the fine detection of cave shape. In this paper, a method of submarine karst morphology detection based on multi-frequency ultrasound is proposed on the basis of laying the ultrasonic probe down the borehole to the rock detection area by cable. Firstly, according to the morphological characteristics of submarine karst, the acoustic reflection model of submarine karst surface is established. Combining with the mechanism of ultrasonic acoustic propagation, the principle of automatic frequency matching is established to ensure that the detection system is always in the best working state, so as to obtain more effective and accurate acoustic detection data. Then, the principle of synthetic aperture imaging is used to realize the various submarine karst caverns. High-precision image of deep horizontal section is presented, and an improved image segmentation method is proposed to extract the contour features of horizontal section at different depths of submarine karst. Search criteria for retrieving calibration value of contour line size are established, and a precise determination system of ultrasonic scanning acoustic beam ranging value is formed to realize the radial distance inversion of the contour of submarine karst horizontal section. Finally, the wheel of submarine karst cavity section is completed. On the basis of fine profile measurement, the depth information of submarine karst is superimposed, and the three-dimensional shape reconstruction of submarine karst is realized by synchronous triangular network method. The feasibility and accuracy of the method are verified by field verification and result analysis. Ó 2020 Elsevier Ltd. All rights reserved.

1. Introduction With the vigorous development of underground space, more and more karst tunnels are constructed. As tunnels are located in various complex hydrogeological rock masses, various karst problems have become the most difficult problems to deal with in the construction of tunnels. Karst is a geological phenomena in which water mainly dissolves soluble rocks, supplemented by mechanical actions, such as water erosion, potential erosion and collapse. The karst development is complex and changeable, with many types, varying sizes and shapes. The main forms of karst development are karst caves, karst troughs, depressions, funnels,

⇑ Corresponding author. E-mail address: [email protected] (X.-h. Huang). https://doi.org/10.1016/j.measurement.2020.107532 0263-2241/Ó 2020 Elsevier Ltd. All rights reserved.

falling caves, collapses and underground rivers [1–4]. The unique geological disaster in karst area is karst collapse, which will have a serious impact on the normal construction and operation of traffic network, and pose a greater threat to the safety of highway and railway. In tunnels constructed in karst areas, mud burst, water burst and landslides often occur, resulting in casualties and serious damage to construction equipment, affecting the construction period and increasing the cost [5,6]. Before the design and construction of the tunnel, the engineering geological conditions of the tunnel engineering section should be investigated in detail. However, due to the inaccuracy of geological exploration and the complexity of karst tunnel geology, the actual surrounding rock conditions revealed after the tunnel excavation may be different from the geological exploration data, which directly leads to the unclear understanding of karst distribution and morphology in

2

J.-c. Wang et al. / Measurement 155 (2020) 107532

the actual construction process. Unexpected accidents such as collapse, roof fall and water gushing often occur in construction, which not only affect the construction period, but also cause damage to construction machinery and casualties [7–10]. Therefore, in order to avoid and deal with all kinds of karst problems encountered, the key is to accurately determine the karst morphology before construction. Submarine karst has a more prominent impact on the construction and design of submarine tunnels. Because of the special fluid structure of the ocean, higher engineering requirements are put forward for the construction of submarine tunnels. Once a submarine tunnel leaks or subsides, its harmful consequences will be incalculable. Therefore, more attention should be paid to the existence of karst in the submarine, and more detailed investigation and detection should be done [11,12]. Because of the concealment, unpredictability and complexity of submarine karst, the commonly used karst detection methods at present, such as seismic method, electromagnetic method and electrical method, have the characteristics of multi-solution and low precision when understanding the size and shape of tunnel karst, which often leads to inaccurate karst detection, which will affect the design and construction of submarine karst to a great extent [13,14]. Improving the accuracy of karst morphology detection has become an urgent problem to be solved in submarine tunnel engineering. Because submarine karst is usually filled with seawater, and the propagation characteristics of ultrasonic wave in seawater are far superior to electromagnetic wave and light wave, we use ultrasonic technology to detect the morphology of submarine karst. Conventional ultrasonic topography mapping technology often shows poor detection effect when detecting objects with irregular reflective surface, and is prone to large errors or problems of nonreflective echo signals. However, submarine karst is formed in natural environment whose wall surface is irregular. It is difficult to obtain high detection accuracy by using conventional ultrasonic topography mapping technology. In view of the shortcomings of traditional ultrasonic measurement technology and the difficulty of fine detection of submarine karst morphology, this paper proposes a method of fine detection of submarine karst morphology based on multi-frequency ultrasound. The principle of automatic frequency matching ensures that the detection system is always in the best working state, so as to obtain more effective and accurate acoustic detection data, and combines with digital images. Processing technology effectively integrates acoustic imaging technology with ranging technology, and realizes fine detection of submarine karst morphology, which can provide more accurate and effective basic data for construction control and safety and stability evaluation of submarine karst, and ensure the construction safety and structural stability of submarine tunnel.

model [15,16]. However, due to the irregular surface of karst body, it is necessary to establish a karst acoustic reflection model considering irregular surface. The acoustic impedance of the transducer is Z1, the acoustic impedance of the medium used as the intermediate layer is Z2, the distance of the medium is h, and the acoustic characteristic impedance of the karst wall is Z3. The transmitted sound wave Pa of the ultrasonic transducer is incident vertically, and the transmitted wave to the transmitting medium is Pb. The transmission wave attenuates to Pc in the propagation medium, and the reflected wave is Pd and the transmitted wave Pt in the rock wall is obtained when incident to the distance h. With the scattering linearity, the reflected wave Pd attenuates to Pe in the propagation medium, and the transmitted wave transmitted by the reflected wave Pe through the surface of the ultrasonic transducer is Pf. Its acoustic reflection model is shown in Fig. 1. The rough surface of rock wall inside submarine karst will cause ultrasonic scattering. On the surface of rock wall detected by local ultrasound, the local area can be approximately regarded as a problem of acoustic scattering with large scale irregular surface. The sound field is calculated by Helmholtz integral. The surface of rock wall detected by ultrasonic wave has rough fluctuation and smoother reflector [17]. Therefore, the acoustic field in the detection area mainly includes two parts: the scattering field generated by rough fluctuation and the effective reflection field generated by smoother reflector, in which the scattering energy will be superimposed on the echo to form the incoherent component. On the surface of the rock wall with smoother reflector, the coherent component of the reflection signal will be formed in the most effective area of the reflection signal. In the case of considering the vertical incidence and reception of ultrasound, the RMS sound pressure value of the received acoustic signal can be expressed as follows [17]:

 1=2 ZZ
1=2 1 Pb 2 rcoh 2 100:2ah þ Pb 2 P2 ¼ ms ðh; /Þ100:2wðh;h;/Þ ds 2h s ð1Þ In the above formula, P is the Rayleigh parameter; h is the distance between the transducer and the rock wall; rcoh is a coherent reflection coefficient and a is a attenuation coefficient; ms (h, u) is a surface scattering coefficient; h is the incident angle of the wave; u is the plane wave parameter. According to the received sound intensity relationship, Ie is the received sound intensity, Ib is the emitted sound intensity, re is the effective reflection coefficient, it can be obtained that:

Ie ¼


r 2 e

2h

Ib eah

ð2Þ

2.2. Optimum working frequency 2. Multifrequency ultrasound detection technology

It can be seen that the karst internal environment from the karst acoustic reflection model, such as the sound velocity of the

2.1. Karst acoustic reflection model When using ultrasonic to detect or measure, the ultrasonic wave needs to be incident into the propagation model composed of different propagation media, and the reflected or transmitted wave obtained by signal processing and analysis can realize the purpose of detection or measurement. When the ultrasonic wave is incident into two different media with interface, the refraction and reflection of the ultrasonic beam will occur at the interface. Among them, one part of the acoustic wave will penetrate the interface and enter another medium to form a transmission wave; the other part of the acoustic wave will be reflected by the interface. The simplest model of plane acoustic wave reflection and transmission in infinite medium is the vertical incident interface

Fig. 1. Acoustic Reflection Model of Submarine Karst.

J.-c. Wang et al. / Measurement 155 (2020) 107532

propagation medium, the size of insoluble particles, the surface fluctuation of the detection object, the reflection intensity of the detection object and the detection distance, which will affect the acoustic propagation characteristics. In different detection environments, in order to simultaneously satisfy the detection range and ensure the detection accuracy, it is necessary to select the optimum working frequency for acoustic detection of karst. In this paper, when using sonar measurement technology to detect karst morphology, the influence of detection distance on the choice of ultrasonic frequency is mainly considered. According to sonar equation, the optimal frequency function of detection distance f ðr Þ is established. If the distance of the detected object remains unchanged, the optimal frequency of the ultrasonic transmitting/ receiving transducer FMcan be deduced according to the quality factor [18], which is defined as:

FM ¼ SL4 =ðNL  DI þ DT Þ

ð3Þ

In the above formula, SL is the sound source level of the ultrasonic probe, NL is the noise level of the surrounding environment of the ultrasonic transducer, DIis the receiving directivity index of the ultrasonic transducer, DTis the detection threshold of the ultrasonic transducer. The specific form of the optimal frequency depends on the dðFMÞ=df calculation results. In order to ensure a fixed beam angle, the receiving directivity index is taken as a constant, which is independent of frequency. It is assumed that the receiving bandwidth is proportional to the working frequency, i.e.Df ¼ Kf , the radiation power is constant, and the detection threshold is independent of the frequency.

dðSLÞ dðDIÞ dðDT Þ ¼ ¼ ¼0 df df df

ð4Þ

When the frequency changes one octave and the change rate of noise level NL is 3dB, taking Df as a constant, the optimal frequency function f ðr Þ can be deduced as:

f ðr Þ ¼

11:6 r 2=3

ð5Þ

According to the relationship (5), when the operating distance r is 1 m~10 m, the optimal frequency f range is 1.2 MHz~250 kHz. In the morphological detection of submarine karst, it is difficult to accurately give the operating distance r between the ultrasonic transducer and the target because the detection scale is unknown, but in order to make the working frequency of the ultrasonic transducer near the optimal detection frequency, the estimated operat

ing distance r can be used. The optimum frequency of the ultrasonic transducer can be calculated separately. The frequency of the transducer varies by one octave and the change rate of noise level is 3 dB. Then the optimum working frequency function f of the transducer can be expressed as follows:


 11:6 f r ¼  r 2=3

ð6Þ 

In the above formula, r is the predicted detection distance between the karst wall and the ultrasonic transducer, which can be obtained by measuring geological data or sonar equipment. 2.3. Principle of multi-frequency ultrasound for submarine karst detection The multi-frequency ultrasonic karst detection system is mainly composed of a borehole probe and a ground system, as shown in Fig. 2(a). The probe in the hole includes cable, sound velocity measurement system, rotary drive mechanism, electronic compass, temperature and pressure sensor, etc. The ground system mainly

3

includes computer and system software, generator, cable car, depth encoder, etc. The working flow of the detection system is as follows: firstly, the multi-frequency ultrasonic probe enters the karst along the drilling wellbore. After the probe is lowered to the specified depth in the submarine karst, the multifrequency ultrasonic probe begins to detect the karst level. After receiving the rotation instruction of the ground control system, the ultrasonic transducer driven by rotary drive mechanism begins to rotate stepwise from the North direction, and rotates at every fixed angle. Horizontal rotation is carried out clockwise or counterclockwise. After each rotation reaches a fixed angle, the ultrasonic transmitter transmits pulse signal to the rock wall, carries out an ultrasonic scanning detection, and the ultrasonic receiver receives the echo signal, and transmits the parameters of the transmitting signal, the echo signal, the azimuth signal and the control signal to the ground system through the cable, and ultrasonic transducer rotates a step angle. The machine rotates a step-by-step angle and continues the above process until the 360° scanning covering the whole horizontal section outline completes the horizontal section scanning of a certain depth of karst. Then, by changing the depth of the probe in the hole, the horizontal section scanning of submarine karst at different depths can be realized. Finally, through data analysis and data processing, the fine exploration of submarine karst morphology can be realized. The working schematic diagram is shown in Fig. 2(b). Multi-frequency ultrasonic detection of submarine karst is a method combining ultrasonic imaging with ultrasonic ranging. Because ultrasonic horizontal section imaging can present the horizontal section characteristics of submarine karst more intuitively and intuitively, as shown in Fig. 3(a), which is the ultrasonic image of a cave horizontal section obtained by the multi-frequency ultrasonic karst detection system. Ultrasonic horizontal section ranging can more accurately reflect and measure the distance length of horizontal scanning line of submarine karst, such as Fig. 3(b) shows that the data fusion of submarine karst horizontal section imaging and ranging is realized by specific analysis algorithm. As shown in Fig. 3(c), the fine inversion of submarine karst horizontal section profile is completed, and the fine detection of submarine karst morphology is realized by combining the multi-level section profile measurement at different depths. The basic principle of submarine karst horizontal section imaging is: when the ultrasonic transducer is located in the submarine karst, the ultrasonic transducer emits acoustic pulse vertically to the karst wall, and uses the ultrasonic receiving transducer to receive the reflected echo from the karst wall, and uses the intensity of the echo to image the horizontal section of the karst, so as to obtain the acoustic impedance information of the rock wall ant the propagation medium in the karst horizontal section. After receiving the echo signal of the horizontal karst section rock wall and propagation medium, the horizontal section imaging of the karst seabed is realized by interpreting the echo signal. Suppose that the ultrasonic wave is projected vertically at a point on the rock surface of horizontal section, and the distance between the point and the ultrasonic transducer is D. The energy reflected from the scanning point by the ultrasonic transducer is F(D), and the output voltage Ub can be expressed as:

U b ¼ 2pI0 KRf expð2aDÞF ðDÞ

ð7Þ

In the formula: I0 is the sound intensity on the surface of the ultrasonic transducer; K is the acoustic-electric conversion coefficient of the ultrasonic transducer; Rf is the acoustic reflection coefficient of the horizontal section reflector of the karst seabed; and a is the propagation loss coefficient of the ultrasonic wave in the propagation medium. The basic principle of submarine karst horizontal section ranging is: when the ultrasonic transducer is located in the submarine

4

J.-c. Wang et al. / Measurement 155 (2020) 107532

(a) Physical Map of Detection System

(b) Working Flow Diagram

Fig. 2. Physical and working diagram of the detection system.

(a) Ultrasound Imaging

(b) Ultrasound Ranging

(c) Profile

Fig. 3. A sketch map of horizontal karst profile detection.

karst, the ultrasonic transducer emits acoustic pulse vertically to the karst wall, and uses the ultrasonic receiving transducer to receive the reflected echo from the karst wall, and calculates the horizontal distance between the ultrasonic transducer and the karst wall according to the time of signal transmitting and receiving. At the same time, the horizontal rotation of the ultrasonic transducer is driven by a stepping motor to realize horizontal cross-section scanning of karst at a certain depth. The sound velocity V of the medium propagating at the known depth of karst at the seabed is measured and the time t of the ultrasonic pulse to and fro is measured. The horizontal distance D between the ultrasonic transducer and the karst wall can be obtained. Its expression is as follows:

D¼

vt 2

ð8Þ

there are a lot of clutter or irrelevant interference signals in the echo signal, the echo signal needs to be processed, while the time-domain correlation processing party needs to process the echo signal. The method is one of the simplest and most effective processing methods in ultrasound imaging [19–21]. For this reason, the time domain correlation processing algorithm is used to process the echo signal. The acoustic signal varies with the distance between the target point of synthetic aperture and the equivalent transducer array element. Usually, this part of calculation is completed in advance, and the results are stored and processed. The working sketch is shown in Fig. 4. In the ultrasonic imaging signal processing of submarine karst horizontal section, there is a good correlation between reflected echoes from target points of karst wall, but the correlation between echo noise is very weak. Therefore, before linear amplitude superposition processing, the adjacent echo signals are multiplied, and then the original echo signals sði; jÞ and sði; j þ 1Þ

3. Fine reconstruction of submarine karst form 3.1. Profile imaging of horizontal section In order to improve the accuracy of profile imaging, more effective data are needed in the process of seabed karst horizontal section profile imaging. In this paper, we make full use of the full sequence of acoustic pulse signals obtained by the multifrequency ultrasonic karst detection system to realize the horizontal profile imaging. Firstly, the acoustic data of each scanning line on the horizontal section are read, then the two-dimensional image of submarine karst horizontal section is composed of the echo intensity of each sampling point in the echo signal. Because

Fig. 4. Sketch of horizontal section contour imaging.

5

J.-c. Wang et al. / Measurement 155 (2020) 107532

received by the transducer at each adjacent position are correlated, and the relationship is expressed by correlation coefficient ci;j . It is assumed that the reflected echo signal matrix of a target point in a deep section of submarine karst is m rows and n columns, in which m represents the sampling points of echo signal in the scanning direction, n represents the scanning times on a section (circumference), and the angle between two adjacent scanning lines is 2*pi/n. Assuming that the ultrasonic transducer receives the echo signals of column j and column j + 1, the correlation coefficients between the two signals can be expressed as follows:

cov½sði; jÞ; sði; j þ 1Þ

ci;j ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi

ð9Þ

Dðsði; jÞÞ  Dðsði; j þ 1ÞÞ

Among them, i 6 m,j þ 1 6 n, and there are the following relations:

cov½sði; jÞ; sði; j þ 1Þ ¼

m h X

ih i sði; jÞ  sði; jÞ sði; j þ 1Þ  sði; j þ 1Þ

i¼1

ð10Þ Dðsði; jÞÞ ¼

m h X

sði; jÞ  sði; jÞ

i2

ð11Þ

i¼1 m h i2 X sði; j þ 1Þ  sði; j þ 1Þ Dðsði; j þ 1ÞÞ ¼

ð12Þ

i¼1

In relation (10), ci;j can be used to evaluate the correlation between two received echoes. The range of variation is from 1 to + 1. When ci;j =+1 indicates that the echoes of column j and column j + 1 are completely (or inversely) linearly correlated. When ci;j is located in between 1 and +1, they have a certain correlation,

 and the non-linear function v ci;j ¼ ci;j exp ci;j is taken as the

formation to segment the image. First, a threshold of each pixel is determined, then the gray value of each pixel is compared with the determined threshold. According to the comparison results, the pixel is regarded as the foreground or background. Geometric contour of rock wall in horizontal section image of submarine karst is a closed area, and the gray value of contour area is different from that of non-contour area. Therefore, this method separates the contour of rock wall from non-contour area by determining appropriate thresholds, and forms a binary image [22,23]. The calculation formula is as follows:

g ðx; yÞ ¼

1 f ðx; yÞ P T 0

Generally, the method of maximum variance between classes is used for image segmentation. Let the gray scale range of horizontal karst section image of seabed be [0, . . ., L  1], and the pixel with gray scale i is ni. Its occurrence probability is pi = ni/N, and N is the total number of pixels. Its value is:

N¼

L1 X

ð16Þ

ni

i¼0

If T is used as the pixel threshold, the horizontal karst section image can be divided into C1 and C2. C1 is composed of pixels with gray value in [0, T  1] and C2 is composed of pixels with gray value in [T, L  1]. Then the probabilities of region C1 and C2 are respectively:

8 T1 P > > pi > < P1 ¼ i¼0

i¼T

The average gray levels of regions C1 and C2 are:

8 T1 P > 1 > ipi > < l1 ¼ P 1

tiplied and the reconstruction points are weighted when they are overlapped. For the weight function, the larger the absolute value is, the stronger the target points in the reconstructed image will be, while the scattering noise points will be suppressed, the signal-to-noise ratio of the image will be enhanced, and the image quality will be greatly improved. After introducing correlation weighting, the algorithm of ultrasonic synthetic aperture imaging for horizontal karst section of seabed can be expressed as follows:

LP 1 > > > ipi : l2 ¼ P12

f ða; bÞ ¼

m X i¼1

v

Pn1 j¼1

n

ci;j

!( n1
Y

v ci;j



   s li;j ; j þ 1

) ð13Þ

j¼1

In the above formula, there are the following relations:



a ¼ 2pr

b ¼ 2pR

ð14Þ

Among them, r is the rotation radius of the transducer placed in the center of the borehole, a represents the circumference of the scanning circle, that is, the circumference of the circular array, R is the radius of the borehole, and b represents the circumference of the circle of the borehole inner wall. 3.2. Profile extraction In order to extract the morphological features of rock wall contour in the horizontal section image of submarine karst, it is necessary to divide the image into several specific and unique regions, i.e. image segmentation. Threshold segmentation is a commonly used image segmentation technology. It uses gray threshold trans-

ð17Þ

LP 1 > > > pi : P2 ¼

weight generating function.


 The adjacent signals are weighted in v ci;j when they are mul-

ð15Þ

f ðx; yÞ < T

i¼0

ð18Þ

i¼T

The average gray level of the whole submarine karst horizontal section image is as follows:

l¼

L1 X

ipi ¼ P1 l1 þ P2 l2

ð19Þ

i¼0

The total variance of the two regions is











rB ¼ P1 l1  l 2 þ P2 l2  l 2 ¼ P1 P2 l1  l2 2

2

ð20Þ

Let T be selected in the range of [0, L  1], so that the rB 2 maximum T value is the optimal region segmentation threshold. The principle of maximum inter-class variance method is simple and the processing speed is fast, but it can not be well adapted to the actual image. There are some differences between the gray value of rock wall contour and the gray value of non-contour area in the horizontal section image of submarine karst. The gray value of adjacent pixels at the boundary of rock wall contour varies greatly, so gradient can be used to express the discontinuity of rock wall contour edge, thus realizing the sharpening of rock wall contour. The horizontal section image of submarine karst is defined as a two-dimensional function f ðx; yÞ, and its gradient is defined as:



rf ¼ gradðf Þ ¼

Gx Gy



ð21Þ

Among them, Gx and Gy is an approximate first-order partial derivative. Its magnitude is:

6

J.-c. Wang et al. / Measurement 155 (2020) 107532

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi jrf ðx; yÞj ¼ Gx 2 þ Gy 2

ð22Þ

In order to sharpen the contour of submarine karst wall and achieve better image segmentation, the gradient operator and the maximum inter-class variance method are combined for image segmentation, that is, the improved variance method. The specific steps of the method are as follows: firstly, the denoised horizontal karst section image F1 is gradient processed to form the gradient image F2, then F2 is binarized to get the binary image F3, and then the corresponding points of the two-dimensional array of images F1 and F3 are multiplied to get a new image F4. Finally, the image F4 is processed by the method of maximum inter-class variance to get the threshold T, which is the optimum threshold that is suitable for image F1. 3.3. Profile reconstruction After the extraction of section contour of karst wall on the seabed is completed, it is necessary to calculate the size of the contour. Because the ultrasonic beam is conical, its reflection characteristics are easily affected by the reflection surface of karst wall. If the center line of the ultrasonic beam intersects with the reflection surface of the protruded karst wall, the first wave of the reflected echo signal received by the ultrasonic receiving transducer is the protruded karst. If the center line of the ultrasonic beam intersects the reflection surface of the depressed karst wall, the reflected echo signal received by the ultrasonic transducer is the reflection signal of both sides of the depressed karst wall, which can not represent the reflection signal of the depressed karst wall directed by the center line of the ultrasonic beam. Therefore, in order to improve the accuracy of the calculation of the size of the contour line, it is necessary to receive the reflected signal of the depressed karst wall. Accurate ultrasonic scanning beam ranging value is used as the calibration value of contour size inversion. Firstly, a plane rectangular coordinate system, called plane image coordinate system, is established based on the position of the probe in a fixed depth rock wall section of submarine karst. The rock wall contour line and the ultrasonic beam scanning line are depicted synchronously, as shown in Fig. 5(a). The angle of the ultrasonic beam is h, the scanning beam is N, and the rock wall contour line intersects with the ultrasonic beam in area D. The outermost intersection point of D clockwise direction in the intersection area is Da, the intersection point of D central angle in the intersection area is Db, the outermost intersection point of D counterclockwise direction in the intersection area is Dc, and the intersection line of DaDc and ODb or its extension line on both sides of the ultrasonic beam intersects Dd, as shown in Fig. 5(b) and (c). It can be seen from Fig. 5(b) and (c) that if ODb is less than ODd, the first reflected echo signal received by the ultrasonic receiving transducer is the signal reflected from the protruded karst wall. If ODb is larger than ODd, the first reflected echo signal received

(a) Ultrasound Imaging

by the ultrasonic receiving transducer is the signal reflected from both sides of the caved karst wall, which cannot represent the depression pointed by the center line of the ultrasonic beam. The reflection signal of karst wall; if the ODb corresponding to all ultrasonic beams at a fixed depth is the smallest, the first wave of reflected echo signal received by the ultrasonic receiving transducer can represent the reflection signal of the karst wall more effectively, and the resolution is stronger. Therefore, by searching for the ultrasonic beams corresponding to the smallest ODb, the optimal one can be found in the section of a fixed depth rock wall. Scanning line sound beam is used as the calibration value of contour size inversion. The calibration value vh of contour size at depth h is set, and its expression is as follows:



vh ¼ min ODh;i;b



ð23Þ

In the above formula, the detection distance of the first wave of the ultrasonic beam corresponding to the scanning line of the i ultrasonic beam at the depth h is indicated as ODh;i;b . In calculating the ODh;i;b corresponding acoustic wave transfer time, the arrival time of the first wave is calculated according to the characteristics of relatively large energy of the first wave beam. In order to minimize the time of real-time operation of searching the maximum energy interval, the energy center convergence method is adopted. Firstly, at depth h, the ultrasonic beam W corresponding to the scanning line of the first ultrasonic beam is divided into four sub-windows, and the energy of each sub-window is calculated. Then the two sub-windows are combined into three new subwindows. The sub-windows corresponding to the largest energy in the new sub-window is W0 . The iteration steps are repeated. After the iteration is terminated, the sub-windows are moved further on the time axis to find out more. Add the exact maximum energy area. According to the principle of mass center, the time T of arrival of the first wave is:

PN Ai t i T ¼ Pi¼1 N i¼1 Ai

ð24Þ

In the above formula, Ai is the amplitude of ti time in the region of the maximum energy corresponding to the ultrasonic beam. At depth h, assuming that the number of pixels in the x-axis direction is Mh and the number of pixels in the y-axis direction is Nh in the plane image coordinate system, according to the corresponding relationship between the calibration value of contour size and the image scanning line, the nh represents the horizontal distance length of a pixel in the section image of karst depth at depth h. The fh represents the length of longitudinal distance, and there exists the following relationship:



nh ¼ vh cosa=Mh ¼ v h T h;i cosa=M h fh ¼ vh sina=Nh ¼ v h T h;i sina=N h

(b) Scanning of Upheaval Rock Wall

(c) Scanning of Depression Wall

Fig. 5. Ultrasound imaging and wall scanning schematic diagram.

ð25Þ

7

J.-c. Wang et al. / Measurement 155 (2020) 107532

In the above formula, the angle between the scanning line of the image and the direction of the x-axis is expressed as a, and the v h represents the propagation speed of sound in the medium. After calibrating the contour size of the optimal scanning line, the contour size on the section image of submarine karst depth can be retrieved according to the length value represented by the pixel points, thus reconstructing the image and size of the section contour. 3.4. Reconstruction of karst stereomorphology Before reconstructing the three-dimensional morphology of submarine karst, it is necessary to establish a space rectangular coordinate system. The center point of vertical boreholes on the surface is the coordinate origin O. According to the right-hand spiral rule, the origin point is selected to point to the north of the geography as the X-axis, the east of the geography as the Y-axis, and the drilling direction is Z-axis. According to the expression of space rectangular coordinate system, the transducer rotates at a uniform speed at the depth h of the borehole, and the rotation angles between the N distance values rh,i and n scanning points are equal. Assuming that the Ph,i point is the i scanning points on the depth h, the polar coordinates of Ph,i point can be transformed into space coordinates. The transformation formula is as follows:


 8 > X h;i ¼ r h;i cos 2pðni1Þ > > <
 2pði1Þ > Y h;i ¼ r h;i sin n > > : Z h;i ¼ h

ð26Þ

On the same vertical section of karst, that is, the scanning angle is equal, the five adjacent points on different depths are Ph2,i, Ph1, i, Ph, Ph,i, Ph+1,i, Ph+2,i, Its corresponding spatial coordinates are (Xh2,i, Yh2,i, Zh2,i), (Xh1,i, Yh1,i, Zh1,i), (Xh,i, Yh,i, Zh,i), (Xh+1,i, Yh+1,i, Zh+1,i), (Xh+2,i, Yh+2,i, Zh+2,i). The direction of the tangent Qh,i at the scanning point Ph,i is determined by the two adjacent scanning points Ph1,i and Ph+1,i. If the tangent Qh,i of the scanning point Ph,i is parallel to the straight line Q0h;i of the over-scanning point Ph1,i and Ph+1,i and the tangent Qh+1,i of the scanning point Ph+1,i is parallel to the straight line Q0hþ1;i of the over-scanning point Ph,i and Ph+2,i. In the space rectangular coordinate system, the tangent intersection point of the two adjacent points Ph,i and Ph+1,i is calculated, which is the control points H0h;i . As shown in Fig. 6(a), the coordinate expression is as follows:

8 x X Y h;i z Z < X ht;iX h;i ¼ Y yht;iY ¼ Z ht;iZ h;i hþ1;i

:

h1;i

xht;i X hþ1;i X hþ2;i X h;i

hþ1;i

¼

h1;i

yht;i Y hþ1;i Y hþ2;i Y h;i

hþ1;i

¼

ð27Þ

h1;i

zht;i Z hþ1;i Z hþ2;i Z h;i

  The coordinate point xht;i ; yht;i ; zht;i is the space rectangular coordinates of control points Hh,i. Similarly, the space rectangular coordinates of control points Hh1,i, Hh+1,i can be solved. If the five adjacent points on the side of a vertical section are Ph2,i+1, Ph1,i+1, Ph,i+1, Ph+1,i+1, Ph+2,i+1, its corresponding spatial coordinates are (Xh2,i+1, Yh2,i+1, Zh2,i+1), (Xh1,i+1, Yh1,i+1, Zh1, i+1), (Xh,i+1, Yh,i+1, Zh,i+1), (Xh+1,i+1, Yh+1,i+1, Zh+1,i+1), (Xh+2,i+1, Yh+2,i+1, Zh+2,i+1). The direction of the tangent Qh,i+1 at the scanning point Ph,i+1 is determined by the two adjacent scanning points Ph1,i+1 and Ph+1,i+1. If the tangent Qh,i+1 of the scanning point Ph,i+1 is parallel to the straight line Q0h;iþ1 of the over-scanning point Ph1,i+1 and Ph+1,i+1 and the tangent Qh+1,i+1 of the scanning point Ph+1,i+1 is parallel to the straight line Q0hþ1;iþ1 of the over-scanning point Ph,i+1 and Ph+2,i+1. In the space rectangular coordinate system, the tangent intersection point of the two adjacent points Ph,i+1 and Ph+1,i+1 is calculated, which is the control points H0h;iþ1 . As shown in Fig. 6(b), the coordinate expression is as follows:

8 x X Y h;iþ1 z Z < X ht;iþ1X h;iþ1 ¼ Y yht;iþ1Y ¼ Z ht;iþ1Z h;iþ1 hþ1;iþ1

:

h1;iþ1

xht;iþ1 X hþ1;iþ1 X hþ2;iþ1 X h;iþ1

hþ1;iþ1

¼

h1;iþ1

yht;iþ1 Y hþ1;iþ1 Y hþ2;iþ1 Y h;iþ1

hþ1;iþ1

¼

h1;iþ1

zht;iþ1 Z hþ1;iþ1 Z hþ2;iþ1 Z h;iþ1

ð28Þ

  The coordinate point xht;iþ1 ; yht;iþ1 ; zht;iþ1 is the space rectangular coordinates of control points Hh,i+1. Similarly, the space rectangular coordinates of control points Hh1,i+1, Hh+1,i+1 can be solved. After the vertical section contour fitting is completed, the points on the horizontal section contour in depth direction are reconstructed according to the triangular mesh model to form the three-dimensional shape of submarine karst. In this paper, a synchronous triangular network model is used to reconstruct the mesh. When connecting points on two adjacent contours, the connection operation is carried out as synchronously as possible on the two contours. It assigns a weight U for each segment of a contour line, which means that the length of the segment is divided by the total length of the contour line where the segment is located. The triangular mesh model is defined as: assuming that UA represents the sum of weights that already exist in the contour line Q hþnh , UB represents the sum of weights that already exist in the contour line Q hþðnh þ1Þ , contour Q hþnh and Q hþðnh þ1Þ are connected to point Phþnh;i ;i and point P hþðnh;i þ1Þ;i respectively, the weights of segment Q hþnh;i ;i and Q hþðnh;i þ1Þ;i are /A and /B respectively, as shown in Fig. 7(a). At this point, the corresponding relationship of the next step is P hþnh;i ;i Phþðnh;iþ1 þ1Þ;iþ1 and Phþðnh;i þ1Þ;i P hþnh;iþ1 ;iþ1 . If

jUA þ /A  UB j < jUB þ /B  UA j, we choose P hþðnh;i þ1Þ;i Phþnh;iþ1 ;iþ1 to move forward along direction A. Otherwise, we choose P hþnh;i ;i Phþðnh;iþ1 þ1Þ;iþ1 to move forward along direction B. According to the above triangular mesh model, we can triangulate the horizontal section contours and the points on the created transition contours in the depth directions of karst exploration to form a more fine karst stereo-morphology map, as shown in Fig. 7(b). 4. Test experiments and analysis

(a) Control point Hh,i

(b) Control point Hh, i+1

Fig. 6. Control Point Solution Diagram.

In order to verify the feasibility and practicability of the fine detection method of submarine karst morphology based on multi-frequency ultrasound, the multi-frequency ultrasound detection system was applied to the measurement of underground karst in a certain sea area of Bohai Sea, China. The whole borehole and submarine karst were scanned in an all-round way, and a large number of detection data were obtained. Through geological drilling and cross-hole CT method, it is estimated that the karst cave is located at the depth of 30–35 m, and its horizontal expansion is about 2 m. The bandwidth of the ultrasonic transducer of the detection system is 100 kHz–1000 kHz.

8

J.-c. Wang et al. / Measurement 155 (2020) 107532

(a) Profile reconstruction

(b) Stereo profile reconstruction

Fig. 7. Schematic diagram of contour shape reconstruction.

4.1. Data acquisition Before the multi-frequency ultrasonic scanning detection is started, the probe of multi-frequency ultrasonic detection system is lowered along the test hole to the position that need to be detected. Then, according to the acoustic reflection model of submarine karst, the detection system chooses the best working frequency to carry out acoustic detection of karst. Because the detection scale is unknown, it is difficult to accurately give the action distance rbetween the ultrasonic transducer and the detection object. But in order to make the working frequency f of the ultrasonic transducer near the optimal detection frequency, the

for scanning beams in different directions of each section is calculated by combining the relation (6). The partial data distribution of the optimal working frequency is shown in Fig. 8. It can be seen from the Fig. 8. The optimal operating frequency is switched to a smaller value in the position with a larger estimated distance, and to a larger value in the position with a smaller estimated distance, so as to ensure that the detection accuracy is higher when the detection range can be reached. By using the optimum working frequency to scan the whole empty area in an all-round way, and upload and save the echo signal, azimuth and depth information of the ultrasonic wave in real time, the efficient acquisition of the detection data of the empty area can be realized.



estimated action distance r can be used to calculate the ultrasound. The optimum frequency of wave transducer is 100 kHz. The distance between the sonar probe and the corresponding karst rock wall is estimated by acquiring 100 kHz acoustic signal. Through the detection of 100 kHz ultrasonic transducer, it is found that there is an empty area between 30.5 and 34.0 m in depth. Seven horizontal sections are scanned in this area, and the depth of each section is estimated. The depth interval of each two adjacent sections is 0.5 m. Some data results are shown in Fig. 8(a). It can be seen from Fig. 8(a) that the size of the estimated distance in the empty area at 30.5–34.0 m varies greatly. In order to switch the working frequency of the ultrasonic transducer to the optimal working frequency state, the optimal working frequency

4.2. Morphological reconstruction After the collection of karst data is completed, the echo signal is processed by time domain correlation processing algorithm, and the contour imaging of each horizontal section of submarine karst is realized by combining the principle of synthetic aperture imaging. The ultrasonic signals of equivalent transducer array elements are stored. When processing the adjacent signals, the weighted processing is carried out, the target points in the reconstructed image are strengthened, and the scattering noise is suppressed. The first horizontal scanning section is processed as shown in Fig. 9(a). An improved variance method combining gradient oper-

(a) Estimating the operating distance

(b) Optimum operating frequency Fig. 8. Part of the probe data graph.

J.-c. Wang et al. / Measurement 155 (2020) 107532

(a) Ultrasound Imaging

9

(b) Profile Profile

Fig. 9. Ultrasound Imaging and Profile.

(a) Contour line intersecting with sound beam

(b) Local enlargement

Fig. 10. Section Profile Curve.

ator and maximum interspecific variance method is used to segment the ultrasonic image and extract the horizontal section contour curve of rock wall, as shown in Fig. 9(b). Based on the section contour curve, the plane image coordinate system is established with the probe location as the origin, and the rock wall contour and the ultrasonic beam scanning line are depicted synchronously, as shown in Fig. 10(a). In order to obtain the optimal calibration value of contour size inversion, the optimal scanning line beam is found by searching the ultrasonic beam corresponding to the minimum ODb. In the first horizontal scanning section, the search result is a scanning line beam with an angle of 100° in the positive East direction. Its local enlargement diagram is shown in Fig. 10(b). The first horizontal scanning section has 575 * 575 pixels, the optimal scanning line has 95 pixels in the x-axis direction, 541 pixels in the y-axis direction and 210 mm in the length of the optimal scanning line. According to the relationship (25), the length of the horizontal distance of the pixels in the first horizontal scanning section is 0.38 mm and the length of the vertical distance is 0.38 mm. It is 0.38 mm, so the actual size of contour curve in each horizontal scanning section can be retrieved. After reconstructing the contour of each horizontal section, the center point of vertical borehole on the surface is chosen as coordinate origin O, the origin point points to the north of geography as X-axis, the east of geography as Y-axis, and the direction of borehole as Z-axis. The spatial rectangular coordinate system is established. The intersection point of ultrasonic scanning beam and the contour of horizontal section is selected as reference modeling

point, and the synchronous triangular network model is adopted, as shown in Fig. 11(a). In order to display the three-dimensional characteristics of submarine karst caves more intuitively and clearly, the optical cave image is matched, as shown in Fig. 11(b). 4.3. Result analysis Comparing the synthetic aperture algorithm with the original signal, it can be seen that the synthetic aperture algorithm can effectively remove the interference noise in the ultrasonic echo signal. Compared with the time-domain correlation processing algorithm, it can be seen that the synthetic aperture algorithm can make the target features in the ultrasonic image clearer. The amplitude of some feature signals is enhanced, which can effectively highlight the profile features of the horizontal section and provide a more detailed image basis for karst cave morphology detection, as shown in Fig. 12. In order to verify the image segmentation effect of the improved variance method used in this paper, it is compared with the maximum variance method and the gradient domain algorithm. From Fig. 13, it can be seen that the maximum variance method is sensitive to the contour information in the image, and the segmentation results are not ideal. The gradient domain algorithm can effectively identify the main contour information, but there are some redundant information. The proposed algorithm will be able to identify the main contour information effectively. The improved variance method, which combines the maximum variance method with the gradient domain algorithm, can sharpen the contours of

10

J.-c. Wang et al. / Measurement 155 (2020) 107532

(a) Non-superimposed image information

(b) Superimposed image information

Fig. 11. Stereo morphological map of submarine karst.

horizontal section images and eliminate the interference of various non-contour information. It also proves that the improved variance method has strong adaptability and can achieve better processing results. In the reconstruction of horizontal section contour, the detection points obtained directly by ultrasonic ranging technology are interpolated and fitted to compare with the image contour extraction method adopted in this paper. In order to show the contrast effect more clearly and intuitively, eight adjacent detection points on the same horizontal section are selected, and the eight detection points are interpolated by polynomial, Spline and linear interpolation. The results of interpolation fitting are shown in Fig. 14. It can be seen from the comparison of interpolation methods that the image contour extraction method can adapt to the situation of different scanning line spacing, while the polynomial interpolation, Spline interpolation and linear interpolation are more suitable for the situation of small and equal scanning line spacing. Therefore, in the reconstruction of horizontal section contour of natural or man-made caves, the image contour extraction method adopted in this paper has stronger applicability, which can provide more accurate and effective reference modeling points for threedimensional modeling of karst, so as to improve the accuracy of karst morphology detection. According to the contour curves of seven horizontal sections and the datum modeling points, the three-dimensional shape grid reconstruction of karst is carried out. As shown in Fig. 15(a), the results show that the depth of karst burial is about 30.5 m, the maximum height of karst is about 3.5 m, the maximum horizontal distance is about 2.7 m, the ratio of height to width of karst is

(a) Original image

about 2.4/2.7 = 0.889, and the shape of karst extends to the Southeast direction. The vertical section diagram of karst extends to the southeast direction as shown in Fig. 15(b), correspondingly. The results of regional cross-hole CT are shown in Fig. 15(c). The corresponding maximum horizontal detection distance is 2.7 m. If only the ultrasonic distance measurement method is used, the maximum detection accuracy is 1%, which means that the maximum resolution of the horizontal section fluctuation based on the single ultrasonic distance measurement method is 27 mm, that is, the response degree of the fluctuation degree is not less than 27 mm. For the fluctuation of the karst contour with the response degree of the fluctuation degree less than 27 mm, the single ultrasonic distance measurement method is difficult to pass the numerical accuracy reflect. The highest resolution of the method described in this paper is 0.38 mm. Theoretically, it can effectively identify and distinguish the features of karst profile with the response degree of fluctuation not less than 0.38 mm. It shows that the method described in this paper has high precision. In combination with the multi-level profile measurement of different depths, it can realize the fine detection of seabed karst morphology. In addition, it can be seen that the basic shape of karst is basically consistent with the karst shape established in this paper from the results of cross-hole CT, but the limit is a cross-section image. It shows that the multi-frequency ultrasonic morphology detection method adopted in this paper can reproduce the threedimensional contour with stronger visuality and more abundant data information. This is mainly due to ultrasound. Wave horizontal section imaging can present the characteristics of karst horizontal section more intuitively and intuitively, while ultrasonic

(b)Temporal correlation processing (c) Synthetic aperture processing Fig. 12. Contrast image of ultrasound.

11

J.-c. Wang et al. / Measurement 155 (2020) 107532

(a) Maximum variance method

(b) Gradient domain algorithm

(c) Improved variance method

Fig. 13. Comparison of the results of contour extraction.

(Red Line-Polynomial Interpolation; Blue Line-Spline Interpolation; Black Line-Linear Interpolation; Dark Blue-Image Contour) Fig. 14. Comparison of contour fitting.

(a) Stereomorphology

(b) Vertical section

(c) Cross-hole CT

Fig. 15. Contrast map of three-dimensional morphology of submarine karst.

horizontal section ranging can more accurately reflect and measure the distance length of horizontal scanning line in karst. By determining the optimal scanning beam, the actual size of contour curve in each horizontal scanning section can be retrieved, and the data of karst horizontal section imaging and ranging can be obtained. By effective fusion, the fine inversion of the horizontal profile of karst is completed, and the fine detection of karst morphology is realized on the basis of the multi-level profile measurement at different depths. 5. Summary Aiming at the difficult technical problems of submarine karst morphology detection, this paper proposes a submarine karst morphology detection method based on the optimum working frequency detection, which integrates the ultrasonic imaging and ranging technology. Firstly, according to the morphological characteristics of submarine karst, the acoustic reflection model of sub-

marine karst surface is established, and the optimal frequency function of the detection system is established based on the mechanism of ultrasonic acoustic propagation, which provides the basis for automatic matching of the ultrasound working frequency. Then, the high-precision maps of submarine karst cavity horizontal sections at different depths are realized by using synthetic aperture imaging principle. An improved image segmentation method is proposed to extract the contour features of submarine karst horizontal sections at different depths. A search criterion for the calibration value of contour line size inversion is established to form a precise determination system of the ultrasonic scanning acoustic beam ranging value, so as to realize the radial distance inversion of the submarine karst horizontal section contour. Finally, on the basis of the precise measurement of the submarine karst cavity section contour, the location information of submarine karst caverns is superimposed, and the three-dimensional shape reconstruction of submarine karst caverns is realized by synchronous triangular network method. The feasibility and accuracy of this

12

J.-c. Wang et al. / Measurement 155 (2020) 107532

method are verified by field verification and result analysis. The results show that: (1) multi-frequency ultrasonic scanning detection technology is feasible and accurate for the shape detection of submarine karst cavity; (2) the optimum working frequency can improve the effectiveness of detection data; (3) synthetic aperture technology can effectively improve the imaging accuracy of ultrasonic horizontal section scanning; (4) the detection accuracy of ultrasonic scanning beam on different rock wall characteristic reflectors exists differences. CRediT authorship contribution statement Jin-chao Wang: Writing - review & editing. Chuan-ying Wang: Data curation. Zeng-qiang Han: Data curation. Yi-teng Wang: Data curation. Xiao-hua Huang: Validation. Declaration of Competing Interest The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. Acknowledgements This work was supported by the State Key Program of National Natural Science of China (Grant No. 41731284); the National Natural Science Foundation for the Youth of China (Grant No. 41902294); the Systematic Project of Guangxi Key Laboratory of Disaster Prevention and Engineering Safety (Grant No. 2019ZDK053). References [1] Wang Muqun, The influence of karst on tunnel engineering and Research on Karst Treatment Technology, Central South University, Hunan, 2011, Master’s degree. [2] Chris Grove, Joe Meiman, Joel Despain, et al., Karst aquifers as at atmospheric carbonsinks: an evolving global network of research sites, us geological survey karst interest group proceedings, Shpherdstown (2002) 20–22.

[3] Pascal Fenart, N.N. Cat, Claude Drogue, et al., Influence of tectonics and neotectonics on the morphogenesis of the peak karst of Halong Bay. Vietnam, Geodinamica Acta (Paris) (12) (1999) 3–4. [4] F. Cabrovsck, W. Dreybrodt, A model of the early evolution of karst aquifers in limestone in the dimensions of lenghth and depth, J. Hydrol. 240 (2001) 206– 224. [5] Wu. Yuehua, Zheng Jixin, Discussion on tunnel karst and treatment technology, Western Explor. Project 3 (9) (2009) 154–155. [6] Mao Yefeng, Wujin, Elementary analysis of control factors and development law of karst development, Western Prospect. Project, 4 (10) (2009) 80–85. [7] Xiong Xuejun, Dayaoshan tunnel disease research and treatment of china railway, China Railway 10 (3) (1997) 39–42. [8] Liang Minggui, Water inflow control of dayaoshan tunnel, World Tunnel 6 (10) (1998) 50–52. [9] Wang Runfu, Sun Guoqing, Li Zhiguo, Construction technology of filling karst cave grouting at the entrance of Yuanliangshan tunnel, Tunnel construction 23 (3) (2003) 28–120. [10] Zhuangjinbo, Water treatment technology of karst pipeline group in Yuanliangshan tunnel, Tunnel Constru. 24 (3) (2004) 20–24. [11] Zhang Dingli, Deformation control technology of unfavorable geological body and structural interface of submarine tunnel, J. Rock Mech. Eng. 26 (11) (2007) 2161–2169. [12] E. Grov, O.T. Blindheim, Risks in adjustable fixed price contracts, Tunnels Tunnell. Int. 35 (6) (2003) 45–47. [13] Wang Jinchao, Wang Chuan Ying, Tang Xinjian, Method for determining the spatial contour of an empty area based on in-hole ultrasound scanning technology, J. Yangtze Acad. Sci. 35 (3) (2018) 104–109. [14] Wang Jinchao, Research on multi-frequency ultrasound scanning detection technology for borehole empty area, Beijing, University of Chinese Academy of Sciences, Master’s degree, 2015. [15] Li Jun, Calculation of matching parameters and location distribution of oil logging transducers, Normal University, Shaanxi, 2009. [16] Quhaolin, Applied Research of Multiphase Interface Measurement System in Crude Oil Storage Tank, Xi’an Petroleum University, 2014. [17] Hang Ruheng, Xu Lianrui, Applicability of ultrasonic reflection logging, Logg. Tech., (6) (1984) 17–25. [18] Chen Hejuan, Ma Baohua, Research on ultrasonic transducer for water detection, Instru. Technol. Sensor, (7) (2001) 9–11. [19] Lan Congqing Cailan, Study on synthetic aperture focused ultrasound imaging method, J. Wuhan Univ. Technol. 18 (1) (1995) 84–87. [20] Sun Baoshen, Zhang Fan, Shen Jianzhong, Study on time domain algorithm of synthetic aperture focused acoustic imaging, J. Acoust. 22 (1) (1997) 42–49. [21] Tian Hao, Tao Liang, Yu Shisheng, et al., Appl. Acoust., 13(5) (1994) 15–17. [22] Jin-chao Wang, Chuan-ying Wang, Chang-qi Zhu, Sheng Hu, Borehole imagebased evaluation of coral reef porosity distribution characteristics, Marine Georesources Geotechnol., 35 (7) (2017) 1008–1017, ISSN:1064-119X. [23] Jinchao Wang, Chuanying Wang, Sheng Hu, Qiang Han, Yiteng Wang, Research on extraction method of structural plane parameters of borehole image, Geotech. Mech. 38 (10) (2017) 3074–3080.