A method for estimating ice layer sound speed and thickness using air wave and ice layer flexural wave
By combining time-frequency analysis of air waves and ice bending waves with a single sensor and a particle swarm optimization algorithm, the problem of difficult equipment deployment for monitoring sea ice characteristics in polar environments was solved, and accurate estimation of ice sound velocity and thickness was achieved.
Patent Information
- Application Number
- CN202510286560.6
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2025-03-12
- Publication Date
- 2025-12-05
- Estimated Expiration
- 2045-03-12
AI Technical Summary
Existing technologies for monitoring sea ice characteristics in polar environments require multiple seismic detector arrays, which are difficult to deploy effectively in harsh environments, and there is a lack of miniaturized and intelligent equipment for estimating ice layer sound velocity and thickness.
Using a single sensor, this method combines time-frequency analysis and particle swarm optimization with air wave and ice layer bending wave analysis to calculate the group velocity of ice layer bending waves from the sound velocity of air waves, thereby retrieving the P-wave velocity, S-wave velocity, ice thickness, and seawater sound velocity.
It enables the monitoring of sea ice characteristics using a single sensor in polar environments, avoiding the need for sensor arrays and providing accurate estimates of ice sound velocity and thickness.
Smart Images

Figure CN120103483B_ABST
Abstract
Description
Technical Field
[0001] This invention belongs to the field of polar acoustic measurement and relates to a method for estimating the sound velocity and thickness of ice layers using air waves and ice layer bending waves. Background Technology
[0002] Existing research indicates that the impacts of global climate change are strongest in the Arctic, which has now become the region with the most intense greenhouse effect on Earth. A key characteristic is the accelerated decline in sea ice coverage; both the extent and average thickness of the ice are decreasing at a faster rate than predicted by climate models. Therefore, a more refined description of the physical processes involved in sea ice models is needed, which involves accurately obtaining parameters such as ice thickness, sound speed, salinity, and temperature.
[0003] Sea ice is a seismic elastic waveguide, with its elastic wave field composed of Lamb wave modes and horizontally polarized transverse wave (SH) modes. Lamb waves are decomposed into two groups: antisymmetric modes An and symmetric modes Sn. n = 0 represents the basic mode, and n = 1, 2, 3 represent higher-order modes. Propagation only occurs at frequencies above their cutoff frequency. The antisymmetric mode produces bending motion, while the symmetric mode produces traction-compression motion. However, the solid-liquid interface causes changes in the boundary conditions. The main difference is the additional generation of a Scholte wave, similar to that propagating on a semi-infinite liquid-solid interface, which behaves similarly to the A0 mode (bending wave) in the low-frequency range and exhibits no dispersion effect in the high-frequency range.
[0004] The propagation of seismic waves within sea ice can be used to monitor sea ice properties such as thickness, sound velocity, Young's modulus, and Poisson's ratio. Traditional seismic methods typically require arrays of dozens or more seismic detectors, but the harsh polar environment limits the deployment of large-scale equipment. Utilizing more intelligent, user-friendly, and miniaturized devices to monitor sea ice properties is a promising approach. Summary of the Invention
[0005] In view of the above-mentioned existing research and technology, the technical problem to be solved by the present invention is to provide a method for estimating the sound velocity and thickness of ice layer using air waves and ice layer bending waves based on a single sensor.
[0006] To address the aforementioned technical problems, the present invention provides a method for estimating the sound velocity and thickness of ice layers using air waves and ice layer bending waves, comprising the following steps:
[0007] Step 1: Collect broadband pulse acoustic signals using a sensor fixed to the ice surface. This sensor includes, but is not limited to, instruments and equipment capable of monitoring particle vibrations, such as seismographs, accelerometers, and microphones. Broadband pulse signals can typically be generated by detonating devices such as detonators, firecrackers, or air cannons on the ice surface to excite bending wave signals propagating within the ice. Record the distance between the detonating device and the ice surface sensor, denoted as R.
[0008] Step 2: Perform time-frequency analysis on the broadband pulse signal acquired by the ice surface sensor. Time-frequency analysis methods generally include short-time Fourier transform and wavelet transform. Wavelet transform is recommended here because it has high frequency resolution in the low-frequency band. Obtain the time spectrum of the signal through time-frequency analysis, denoted as P(t,f), where P is the amplitude, t represents time, and f represents frequency. Determine the arrival time of the first arriving curved wave and the subsequent air wave in the time spectrum, denoted as t1.
[0009] Step 3: Measure the temperature of the experimental site and calculate the propagation speed of the air wave, denoted as c1; the relationship between the propagation speed of air and temperature satisfies:
[0010] c1 = 331.6 + 0.6T
[0011] Where T represents temperature, and the unit is degrees Celsius (°C).
[0012] Step 4: Convert the time spectrum to a group velocity spectrum and extract the group velocity dispersion curve of the curved wave; in Step 2, the two coordinates of the time spectrum P(t,f) are time t and frequency f, and the two coordinates of the group velocity spectrum P(c,f) are group velocity c and frequency f. The conversion relationship between the two coordinates of time t and group velocity is as follows:
[0013]
[0014] The correspondence between the time spectrum and the group velocity spectrum is obtained as follows:
[0015] Where c1 is the propagation speed of the air wave, t1 represents the time when the air wave first arrives, R represents the distance between the detonation device and the ice surface sensor, t represents time, f represents frequency, and P(c,f) is the group velocity spectrum.
[0016] Extract the flexural wave from the group velocity spectrum P(c,f) to form a dispersion curve, denoted as c2(f), where f is the frequency point of the flexural wave, and c2(f) represents the group velocity coordinate corresponding to the flexural wave in the group velocity spectrum P(c,f) when the frequency coordinate is f.
[0017] Step 5: Initialize the longitudinal wave velocity c of the ice layer s transverse wave velocity c t Ice thickness d, seawater sound speed cw The initial ranges for these four parameters are not fixed and can be modified using historical measurement data. The default settings are: initial range for P-wave velocity is 3000m / s-4500m / s, initial range for S-wave velocity is 1600m / s-2000m / s, initial range for ice thickness is 0.1m-10.0m, and initial range for seawater velocity is 1400m / s-1500m / s.
[0018] Step 6: Calculate the desired group velocity curve based on the sound propagation model and construct the cost function; the cost function is:
[0019]
[0020] Where N represents the number of frequency points of the dispersion curve extracted in step 4, c2(f n ) is when the nth frequency is f n The group velocity of the time-bending wave. c(f) n The frequency f is calculated using the sound propagation model. n The theoretical group velocity of the curved wave. The acoustic propagation model used here is the Krakel model from Kraken, a normal wave model for underwater acoustic propagation. This model can calculate the horizontal wave number k of the curved wave under given conditions of P-wave velocity, S-wave velocity, ice thickness, and seawater sound velocity at a given frequency f. The relationship between the horizontal wave number and the group velocity is:
[0021]
[0022] in It represents the derivative of frequency f with respect to the horizontal wave number k.
[0023] Step 7: Optimize the cost function using the Particle Swarm Optimization (PSO) algorithm and output the final parameter estimates. PSO is a computationally convenient and fast-converging global optimization algorithm. Based on the parameters to be inverted set in Step 5, the number of particle coordinates is set to 4, and the total number of particles is 100. The upper and lower limits of the coordinates are the initial settings from Step 5. The cost function from Step 6 is minimized using the PSO algorithm, continuously updating the four parameter coordinates of the particles until convergence is achieved. The P-wave velocity, S-wave velocity, ice thickness, and seawater sound speed are estimated based on the average of the 100 particle coordinates.
[0024] The beneficial effects of this invention are as follows: This invention can be implemented using only a single sensor, without the need for an array of multiple sensors; this method does not require the transmitting and receiving devices to be synchronized; it proposes to calculate the group velocity of ice bending waves using the speed of sound and time of arrival of air waves; the relationship between group velocity and acoustic parameters can be used to invert the longitudinal wave velocity, transverse wave velocity, ice thickness, and seawater sound velocity; and its effectiveness was finally verified through experiments on ice. Attached Figure Description
[0025] Figure 1 This is a flowchart of the parameter inversion process in this invention;
[0026] Figure 2 This is a diagram of the field experiment scenario in this invention;
[0027] Figure 3 This is the time spectrum of the signal in this invention;
[0028] Figure 4 This refers to the group velocity spectrum in this invention;
[0029] Figure 5 This refers to the change in the cost function of particle swarm optimization in this invention;
[0030] Figure 6 These are the parameter estimation results in this invention;
[0031] Figure 7 This is the result of the comparison of flexural wave group velocities in this invention. Detailed Implementation
[0032] The present invention will be further described below with reference to the accompanying drawings and embodiments.
[0033] like Figure 1 As shown, the present invention provides a method for estimating the sound velocity and thickness of ice layers using air waves and ice layer bending waves, comprising the following steps:
[0034] like Figure 2 The image shows an outdoor experimental scenario, including a sensor fixed to the ice surface and an explosion sound source that can generate air waves and ice bending waves.
[0035] Step 1 involves acquiring broadband pulsed sound signals using a sensor fixed to the ice surface. This sensor includes, but is not limited to, instruments and equipment capable of monitoring particle vibrations, such as seismographs, accelerometers, and microphones. Broadband pulsed signals are typically generated by detonating devices such as detonators, firecrackers, or air cannons on the ice surface to excite bending wave signals propagating within the ice. The distance between the detonating device and the ice surface sensor is recorded and denoted as R. In the ice experiment, firecrackers were detonated on the ice surface, and a vertical accelerometer was fixed to the ice surface to collect vibration signals at a distance of R = 200 m. The bending wave from the ice surface arrives first, followed by the air wave, which has a larger amplitude than the bending wave.
[0036] like Figure 3As shown, step 2 performs time-frequency analysis on the explosion sound collected by the ice surface sensor. Time-frequency analysis methods generally include short-time Fourier transform and wavelet transform; wavelet transform is used here because it has high frequency resolution in the low-frequency band. The time-frequency spectrum of the signal is obtained through time-frequency analysis, denoted as P(t,f), where P is the amplitude, t represents time, and f represents frequency. The time spectrum shows the arrival of the bending wave first, followed by the air wave, and the arrival time of the air wave is determined, denoted as t1 = 0.665s.
[0037] Step 3: Measure the temperature of the experimental site and calculate the air wave propagation speed, denoted as c1; the relationship between air propagation speed and temperature satisfies:
[0038] c1 = 331.6 + 0.6T
[0039] Where T represents temperature, and the unit is degrees Celsius (°C). The ambient temperature is approximately -20°C, therefore the speed of air travel is approximately 320 m / s.
[0040] like Figure 4 As shown, step 4 converts the time spectrum into a group velocity spectrum and extracts the group velocity dispersion curve of the curved wave; in step 2, the two coordinates of the time spectrum P(t,f) are time t and frequency f, and the two coordinates of the group velocity spectrum P(c,f) are group velocity c and frequency f. The transformation relationship between the two coordinates of time t and group velocity is as follows:
[0041]
[0042] The correspondence between the time spectrum and the group velocity spectrum is obtained as follows:
[0043] Where c1 is the propagation speed of the air wave, t1 represents the time when the air wave first arrives, R represents the distance between the detonation device and the ice surface sensor, t represents time, f represents frequency, and P(c,f) is the group velocity spectrum.
[0044] Extract the curved wave from the group velocity spectrum P(c,f) to form a dispersion curve, denoted as c2(f), as follows. Figure 4 As shown by the dashed line, f is the frequency of the flexural wave, and c2(f) represents the group velocity coordinate of the flexural wave in the group velocity spectrum P(c,f) when the frequency coordinate is f.
[0045] Step 5: Initialize the longitudinal wave velocity c of the ice layer s transverse wave velocity c t Ice thickness d, seawater sound speed c wThe initial ranges for these four parameters are not fixed and can be modified using historical measurement data. The default settings are: initial range for P-wave velocity is 3000m / s-4500m / s, initial range for S-wave velocity is 1600m / s-2000m / s, initial range for ice thickness is 0.1m-10.0m, and initial range for seawater velocity is 1400m / s-1500m / s.
[0046] Step 6: Calculate the desired group velocity curve based on the sound propagation model and construct the cost function; the cost function is:
[0047]
[0048] Where N represents the number of frequency points of the dispersion curve extracted in step 4, c2(f n ) is when the nth frequency is f n The group velocity of the time-bending wave. c(f) n The frequency f is calculated using the sound propagation model. n The theoretical group velocity of the curved wave. The acoustic propagation model used here is the Kraker model from the Kraken standard wave model for underwater acoustic propagation, which can calculate the horizontal wave number k of the curved wave given the P-wave velocity, S-wave velocity, ice thickness, and seawater sound velocity. The relationship between the horizontal wave number and the group velocity is as follows:
[0049]
[0050] in It represents the derivative of frequency f with respect to the horizontal wave number k.
[0051] like Figure 5 As shown, step 7 uses the particle swarm optimization (PSO) algorithm to optimize the cost function and outputs the final parameter estimates. PSO is a computationally convenient and fast-converging global optimization algorithm. Based on the parameters to be inverted set in step 5, the number of particle coordinates is set to 4, and the total number of particles is 100. The upper and lower limits of the coordinates are the initial settings from step 5. The cost function in step 6 is minimized using the PSO algorithm, continuously updating the four parameter coordinates of the particles until convergence is achieved. The probability density of the 100 particle coordinates is used to estimate the P-wave velocity, S-wave velocity, ice thickness, and seawater sound speed. Figure 5 The results show that the optimization function value of particle swarm optimization decreases continuously with the number of iterations, and the final cost function tends to a stable minimum value, indicating that the coordinate parameters of the particles are now consistent with the actual parameters.
[0052] like Figure 6 As shown, the final coordinate parameters of the particles are statistically analyzed, probability densities are plotted, and the relative relationships between the four parameters are displayed. The final estimated parameter values are: the sound velocity of the P-wave in the ice layer, c. s= 4195 m / s, shear wave velocity c t =1712m / s, speed of sound in seawater c w =1434m / s, ice layer thickness d=0.77m.
[0053] like Figure 7 As shown in the figure, the dashed line represents the theoretical group velocity calculated using the parameter values estimated in step 7. A consistency comparison was made between the theoretical group velocity and the observed group velocity based on the experimental observation values given in step 4, showing that the theoretical group velocity is consistent with the observed group velocity, reflecting the accuracy of the parameter estimation method.
[0054] It should be understood that the application of the present invention is not limited to the examples above. Those skilled in the art can make improvements or modifications based on the above description, and all such improvements and modifications should fall within the protection scope of the appended claims.
Claims
1. A method for estimating the sound velocity and thickness of ice layers using air waves and ice layer bending waves, characterized in that, Includes the following steps: Step 1: Use a sensor fixed to the ice surface to collect broadband pulse acoustic signals; record the distance between the detonation device and the ice surface sensor, and denote it as R; Step 2: Perform time-frequency analysis on the broadband pulse signal collected by the ice surface sensor; obtain the time spectrum of the signal through time-frequency analysis, denoted as P(t,f), where P is the amplitude, t represents time, and f represents frequency; Step 3: Measure the temperature of the experimental site and calculate the propagation speed of the air wave, denoted as c1; the relationship between the propagation speed of air and temperature satisfies: c1 = 331.6 + 0.6T Where T represents temperature, and the unit is degrees Celsius (°C); Step 4: Convert the time spectrum to a group velocity spectrum and extract the group velocity dispersion curve of the flexural wave; in Step 2, the two coordinates of the time spectrum P(t,f) are time t and frequency f, and the two coordinates of the group velocity spectrum P(c,f) are group velocity c and frequency f; the transformation relationship between the two coordinates of time t and group velocity is as follows: The correspondence between the time spectrum and the group velocity spectrum is obtained as follows: Where c1 is the propagation speed of the air wave, t1 represents the time when the air wave first arrives, R represents the distance between the detonation device and the ice surface sensor, t represents time, f represents frequency, and P(c,f) is the group velocity spectrum. Extract the flexural wave from the group velocity spectrum P(c,f) to form a dispersion curve, denoted as c2(f), where f is the frequency point of the flexural wave, and c2(f) represents the group velocity coordinate corresponding to the flexural wave in the group velocity spectrum P(c,f) when the frequency coordinate is f. Step 5: Initialize the longitudinal wave velocity c of the ice layer s transverse wave velocity c t Ice thickness d, seawater sound speed c w ; Step 6: Calculate the desired group velocity curve based on the sound propagation model and construct the cost function; Step 7: Optimize the cost function using the particle swarm optimization algorithm and output the final parameter estimates; based on the parameters to be inverted set in Step 5, set the number of particle coordinates to 4 and the number of particles to 100, with the upper and lower limits of the coordinates being the initial settings in Step 5; minimize the cost function of Step 6 using the particle swarm optimization algorithm, continuously updating the coordinates of the four parameters of the particles until convergence is achieved, and estimate the P-wave velocity, S-wave velocity, ice thickness, and seawater sound speed based on the mean of the coordinates of the 100 particles.
2. The method for estimating the sound velocity and thickness of ice layers using air waves and ice layer bending waves according to claim 1, characterized in that, The sensors in step 1 are seismographs, accelerometers, or microphones; broadband pulse signals are detonated on the ice surface using detonators, firecrackers, or air cannons to excite bending wave signals propagating in the ice layer.
3. The method for estimating the sound velocity and thickness of ice layers using air waves and ice layer bending waves according to claim 1, characterized in that, In step 2, the time-frequency analysis method uses wavelet transform.
4. The method for estimating the sound velocity and thickness of ice layers using air waves and ice layer bending waves according to claim 1, characterized in that, In step 5, the initial ranges of the four parameters are modified using historical measurement data; the default settings are: the initial range for P-wave velocity is 3000m / s-4500m / s, the initial range for S-wave velocity is 1600m / s-2000m / s, the initial range for ice thickness is 0.1m-10.0m, and the initial range for seawater velocity is 1400m / s-1500m / s.
5. The method for estimating the sound velocity and thickness of ice layers using air waves and ice layer bending waves according to claim 1, characterized in that, The cost function in step 6 is: Where N represents the number of frequency points of the dispersion curve extracted in step 4, c2(f n ) is when the nth frequency is f n The group velocity of a curved wave; c(f) n The frequency f is calculated using the sound propagation model. n The theoretical group velocity of the curved wave is calculated. The acoustic propagation model uses the Kraker model from the Kraken standard wave model for underwater acoustic propagation. The horizontal wave number k of the curved wave is calculated at a given frequency f, considering the P-wave velocity, S-wave velocity, ice thickness, and seawater acoustic velocity. The relationship between the horizontal wave number and the group velocity is as follows: in It represents the derivative of frequency f with respect to the horizontal wave number k.
Citation Information
Patent Citations
On-ice seismic source positioning method based on bending waves
CN113687308A
Method for radar-location determination of ice thickness
RU2526222C1