CN114152329A - Method for detecting spectrum peak of underwater acoustic signal - Google Patents
Method for detecting spectrum peak of underwater acoustic signal Download PDFInfo
- Publication number
- CN114152329A CN114152329A CN202111411551.3A CN202111411551A CN114152329A CN 114152329 A CN114152329 A CN 114152329A CN 202111411551 A CN202111411551 A CN 202111411551A CN 114152329 A CN114152329 A CN 114152329A
- Authority
- CN
- China
- Prior art keywords
- length
- gaussian
- win
- peak
- waveform
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01H—MEASUREMENT OF MECHANICAL VIBRATIONS OR ULTRASONIC, SONIC OR INFRASONIC WAVES
- G01H11/00—Measuring mechanical vibrations or ultrasonic, sonic or infrasonic waves by detecting changes in electric or magnetic properties
- G01H11/06—Measuring mechanical vibrations or ultrasonic, sonic or infrasonic waves by detecting changes in electric or magnetic properties by electric means
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
- G01R23/00—Arrangements for measuring frequencies; Arrangements for analysing frequency spectra
- G01R23/16—Spectrum analysis; Fourier analysis
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/14—Fourier, Walsh or analogous domain transformations, e.g. Laplace, Hilbert, Karhunen-Loeve, transforms
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/18—Complex mathematical operations for evaluating statistical data, e.g. average values, frequency distributions, probability functions, regression analysis
Abstract
The invention discloses a method for detecting a spectral peak of an underwater sound signal, which comprises the steps of eliminating polarity of a signal waveform containing the spectral peak to obtain B (n), solving a gain coefficient c (n), multiplying the c (n) with the B (n) to obtain an output signal C (n), filtering and smoothing a median of the C (n) to obtain D (n), fitting the D (n) by using a Gaussian window to obtain E (n), and extracting a maximum value point and a corresponding maximum value from the E (n). The method can detect the spectrum peak existing in the signal waveform through the steps of absolute value depolarization, double-sliding-window spectrum peak enhancement, sliding average filtering, Gaussian fitting and the like, effectively reduces the influence of a hydroacoustic channel on the original spectrum peak, and realizes the enhancement of the spectrum peak characteristic.
Description
Technical Field
The invention belongs to the technical field of feature extraction of underwater acoustic signals, and relates to a method for detecting a spectral peak of an underwater acoustic signal.
Background
The detection of a spectral peak is often required when performing underwater acoustic signal identification or parameter estimation, for example, the detection of a spectral peak of a frequency spectrum is required when performing underwater acoustic Frequency Shift Keying (FSK) identification, the detection of an autocorrelation spectral peak is required when performing underwater acoustic Orthogonal Frequency Division Multiplexing (OFDM) signal identification, and the detection of a spectral peak position of a frequency domain is required when performing carrier frequency estimation on communication signals such as FSK, Phase Shift Keying (PSK), Direct Sequence Spread Spectrum (DSSS), and the like.
The spectral peak detection method of signals in underwater sound mostly originates from radio. There are many classical spectral peak extraction methods in radio, such as direct spectral peak detection algorithm, three-point spectral peak detection algorithm, gaussian fitting algorithm, polynomial fitting algorithm, etc. Direct spectral peak detection has poor noise immunity, although computational complexity is low. The three-point spectral peak detection, Gaussian fitting and polynomial fitting algorithms have higher requirements on the shape of the spectral peak. The hydro-acoustic channel is more complex than a common wireless channel, has a complex time-space-frequency variation characteristic, shows frequency selective fading in a frequency domain, and a spectral peak obtained by using correlation also becomes unstable when a signal-to-noise ratio is low or a channel time is significantly changed. Spectral peak detection using the signal frequency domain or correlation waveform in the hydroacoustic channel is highly challenging.
Disclosure of Invention
In view of the foregoing prior art, the technical problem to be solved by the present invention is to provide a method for detecting a spectral peak of an underwater acoustic signal with an unobvious spectral peak feature in an underwater acoustic channel, which can effectively enhance the spectral peak feature and suppress interference.
In order to solve the technical problem, the method for detecting the spectrum peak of the underwater sound signal comprises the following steps:
s1, a (N) represents a signal waveform including a spectral peak, N is the position of the nth sampling point of the signal waveform, N is 1,2, N is the number of the sampling points, when the spectral peak waveform has negative polarity, b (N) is | a (N) |, when the spectral peak waveform has no negative polarity, b (N) is a (N) |;
s2, establishing two sliding windows with the nth sampling point of B (n) as the center, wherein the length of the short sliding window is 2Swin+1, long sliding window length 2lwin+1, and lwin＝2swin(ii) a Obtaining a gain factor c (through the energy difference of the two sliding windows)n)，
Is the short sliding window energy of the nth sample point,
is the long sliding window energy of the nth sample point, where
The output signal is C (N) ═ B (N) · c (N), where N ═ 1,2, ·, N, when c (N) is negative, N-swin≤0、n-lwin≤0、n+swin> N or N + lwinIf > N, let c (N) be 0;
s3, performing moving average filtering on the c (n) waveform obtained in S2, and obtaining an output signal after moving average filtering:
wherein 2P +1 is the moving average filter length, N ═ 1,2, ·, N, when N-P ≦ 0 or N + P > N, let d (N) ═ c (N);
s4, performing Gaussian fitting on the D (n) obtained in S3, wherein the Gaussian window function is as follows:
wherein L is the length of the Gaussian window function, - (L-1)/2 ≦ n ≦ L-1)/2, α is the width coefficient of the Gaussian window function, and α is inversely proportional to the standard deviation σ of the Gaussian function (L-1)/2 α; the output signal after gaussian fitting is:
wherein, L ═ 2Q +1 is the length of Gaussian window function, N ═ 1,2, ·, N, then normalize E (N) to get the optimized spectrum peak waveform, when N-Q ≤ 0 or N + Q > N, let E (N) ═ D (N);
s5, extracting all the maximum values and the corresponding maximum value points in the E (n), and sorting the maximum values from large to small, if the spectrum peak is a single spectrum peak, extracting the maximum value and the corresponding maximum value point in the E (n) as the spectrum peak, if the spectrum peak has a plurality of spectrum peaks, extracting the top m maximum values and the corresponding maximum value points in the E (n) as the spectrum peaks, wherein m is equal to the number of the spectrum peaks.
Further, the moving average filter length in S3 is equal to the short sliding window length in S2.
Further, the length of the gaussian window function in S4 is equal to the length of the short sliding window in S2.
Further, the length of the short sliding window in S2, the length of the moving average filter in S3, and the length of the gaussian window function in S4 are equal.
The invention has the beneficial effects that: the method can detect the spectrum peak existing in the signal waveform through the steps of absolute value depolarization, double-sliding-window spectrum peak enhancement, sliding average filtering, Gaussian fitting and the like, effectively reduces the influence of a hydroacoustic channel on the original spectrum peak, and realizes the enhancement of the spectrum peak characteristic. The invention has two main applications:
in a first aspect, the invention can be used for extracting spectral peak features of non-cooperative underwater acoustic communication signals, and the identification of the non-cooperative underwater acoustic communication signals is realized through the features. For example, in underwater sound OFDM, a cyclic prefix structure generally exists, and we perform spectral peak detection on an autocorrelation waveform of a signal to determine whether two correlation peaks exist. For another example, the underwater sound FSK signal has a plurality of spectral peaks in a frequency domain, spectral peak detection is carried out on a signal frequency domain waveform to determine the number of spectral peaks, and the modulation order of FSK is determined according to the number of detected spectral peaks.
In a second aspect, the present invention can be used to estimate the carrier frequencies of an underwater acoustic FSK signal. The spectrum peak detection is carried out by using the method, and finally, a maximum point and a corresponding maximum value are obtained, and the frequency of the signal is corresponding to the maximum point on a frequency domain.
Drawings
FIG. 1 is a flow chart of a method for extracting spectral peaks of an underwater acoustic signal according to the present invention;
FIG. 2 is a diagram of channel impulse response for testing provided by the present invention;
FIG. 3 is an original frequency domain waveform obtained by Fourier transform of a simulated 8 FSK;
FIG. 4 is a graph of gain coefficients obtained based on the original frequency domain waveform of FIG. 3;
FIG. 5 is a result of multiplying the gain factor of FIG. 4 with the original frequency domain waveform of FIG. 3;
FIG. 6 shows the result of moving average filtering of the waveform of FIG. 5;
FIG. 7 is a result of a Gaussian fit to the waveform of FIG. 6;
fig. 8 is a result of defining the abscissa in fig. 7 within 10kHz to 14 kHz.
Detailed Description
The invention is further described with reference to the drawings and the specific embodiments in the following description.
With reference to fig. 1, the present invention comprises the following steps:
s1, a (n) represents a signal waveform including a spectral peak, where n is the position of the nth sample point of the signal waveform. In order to eliminate uncertainty of bipolar signals, an absolute value of a (n) is obtained to obtain b (n), that is, b (n) ═ a (n) |, and when no negative polarity exists in a spectrum peak waveform, b (n) ═ a (n) is obtained;
s2, performing the following operations for all the sampling points one by one: taking the sampling point n of B (n) as the center, establishing a long sliding window and a short sliding window, wherein the length of the short sliding window is 2swin+1, the length of the long sliding window is 2lwin+1, and lwin＝2swin. The gain coefficient c (n) is obtained by the energy difference of the two sliding windows, and can be expressed as
Wherein
When c (n) is negative, n-swin≤0，n-lwin≤0，n+swin> N or N + lwinIf > N, c (N) is 0, and the output signal is c (N) b (N) c (N). By multiplying a gain coefficient c (n), when a spectrum peak exists in the signal waveform, the waveform at the spectrum peak is amplified, and the waveforms at other places are suppressed;
s3, performing the following operations for all the sampling points one by one: the c (n) waveform obtained in S2 is subjected to moving average filtering, which can reduce random noise and approximately obtain the envelope shape of the spectral peak waveform. The output signal after the moving average filtering can be expressed as
Wherein 2P +1 is the length of the moving average filter, when N-P is less than or equal to 0 or N + P is more than N, let D (N) become C (N);
s4, performing the following operations for all the sampling points one by one: and finally, performing Gaussian fitting on the D (n) obtained in the step S3, wherein the step can effectively smooth the D (n), and is favorable for extraction of the maximum values in the D (n) and detection of spectral peaks. The Gaussian window function can be expressed as
Where L is the length of the gaussian window function, - (L-1)/2 ≦ n ≦ L-1)/2, α is the width coefficient of the gaussian window function, and α is inversely proportional to the standard deviation σ of the gaussian function (L-1)/2 α. The output signal after Gaussian fitting can be expressed as
Wherein, L2Q +1 is the length of the Gaussian window function, when N-Q is less than or equal to 0 or N + Q is more than N, E (N) D (N) is made, and then E (N) is normalized to obtain the optimized spectrum peak waveform;
s5, extracting all the maximum values and the corresponding maximum value points in E (n), and sorting the maximum values from large to small according to the sizes of the maximum values. If the spectrum peak is a single spectrum peak, extracting the maximum value and the corresponding maximum value point in E (n) as the spectrum peak, and if a plurality of spectrum peaks exist, extracting the first m maximum values and the corresponding maximum value points in E (n) as the spectrum peaks, wherein m is equal to the number of the spectrum peaks.
The spectral peak refers to a frequency domain spectral peak or a signal correlation spectral peak, the bipolar signal in S1 refers to a case where a signal waveform has a negative polarity due to a correlation operation or the like, and in some cases, there is no case where a negative polarity exists, for example, the signal frequency domain waveform, and the length of the sliding window should be as long as possible to include a spectral peak in the setting of the short sliding window length in S2, but should not be too long, and when a plurality of spectral peaks exist in a signal waveform, the sliding window length is too long, and a plurality of spectral peaks are likely to be included in the window. In addition, the interference suppression capability is reduced when the sliding window length is too short. The actual application can be set according to empirical values. In S3, a large number of pseudo envelopes are easily formed when the length of the moving average filter is too short, and a large amount of useless information is included when the length is too long, and the length of the moving average filter can be set to be consistent with the length of the short sliding window in S2 in order to ensure that the envelope of the spectral peak can be effectively extracted. The length of the gaussian window function in S4 can be kept consistent with the short sliding window in S2 and the moving average filter length in S3.
Examples are given below with specific parameters:
with reference to fig. 1, it is assumed that the acquired signal is an 8FSK signal, the sampling rate is 48kHz, the bandwidth of the 8FSK signal is 4kHz, the carrier frequencies are {10.25kHz, 10.75kHz, 11.25kHz, 11.75kHz, 12.25kHz, 12.75kHz, 13.25kHz, 13.75kHz }, the simulated signal-to-noise ratio is set to be 0dB over the full band, the channel impulse response used in the simulation is shown in fig. 2, and we first perform fourier transform on the acquired 8FSK signal to obtain a signal waveform a (n) including spectral peaks, and the signal waveform of a (n) is shown in fig. 3. The specific implementation of the invention comprises the following steps:
s1, since there is no negative polarity in a (n), the absolute value of depolarization may be omitted, i.e. b (n) ═ a (n);
s2, setting the length of the short sliding window as the length of the sampling point corresponding to the 200Hz bandwidth, respectively calculating the energy in the two sliding windows, and obtaining the gain coefficient c (n) according to the energy difference. The waveform of the gain coefficient c (n) is shown by a broken line in fig. 4. The final output signal is c (n) ═ b (n) · c (n), and the waveform of c (n) is shown in fig. 5;
s3, performing moving average filtering on the C (n) obtained in the S2 to obtain a filtered signal D (n), wherein the waveform of the D (n) is shown in FIG. 6;
s4, performing Gaussian fitting on the D (n) obtained in S3, firstly generating a Gaussian window G (n), fitting each sampling point in D (n) by using G (n) to obtain a fitted signal E (n), and normalizing E (n) to obtain a waveform as shown in FIG. 7.
S5, extracting the maximum points and their corresponding maximum values in fig. 7 as shown in fig. 8, 8 distinct spectral peaks can be found in fig. 8, fig. 8 is a signal waveform in which the abscissa of fig. 7 is limited to the range of 10kHz to 14kHz, and it can be found that the frequencies corresponding to these maximum points can correspond to 8 carrier frequencies of the 8FSK signal.
Claims (4)
1. A method for detecting a spectral peak of an underwater sound signal is characterized by comprising the following steps:
s1, a (N) represents a signal waveform including a spectral peak, N is the position of the nth sampling point of the signal waveform, N is 1,2, N is the number of the sampling points, when the spectral peak waveform has negative polarity, b (N) is | a (N) |, when the spectral peak waveform has no negative polarity, b (N) is a (N) |;
s2, establishing two sliding windows with the nth sampling point of B (n) as the center, wherein the length of the short sliding window is 2Swin+1, long sliding window length 2lwin+1, and lwin＝2swin(ii) a Obtaining a gain coefficient c (n) through the energy difference of the two sliding windows,
is the short sliding window energy of the nth sample point,
is the long sliding window energy of the nth sample point, where
The output signal is C (N) ═ B (N) · c (N), where N ═ 1,2, ·, N, when c (N) is negative, N-swin≤0、n-lwin≤0、n+swin> N or N + lwinIf > N, let c (N) be 0;
s3, performing moving average filtering on the c (n) waveform obtained in S2, and obtaining an output signal after moving average filtering:
wherein 2P +1 is the moving average filter length, N ═ 1,2, ·, N, when N-P ≦ 0 or N + P > N, let d (N) ═ c (N);
s4, performing Gaussian fitting on the D (n) obtained in S3, wherein the Gaussian window function is as follows:
wherein L is the length of the Gaussian window function, - (L-1)/2 ≦ n ≦ L-1)/2, α is the width coefficient of the Gaussian window function, and α is inversely proportional to the standard deviation σ of the Gaussian function (L-1)/2 α; the output signal after gaussian fitting is:
wherein, L ═ 2Q +1 is the length of Gaussian window function, N ═ 1,2, ·, N, then normalize E (N) to get the optimized spectrum peak waveform, when N-Q ≤ 0 or N + Q > N, let E (N) ═ D (N);
s5, extracting all the maximum values and the corresponding maximum value points in the E (n), and sorting the maximum values from large to small, if the spectrum peak is a single spectrum peak, extracting the maximum value and the corresponding maximum value point in the E (n) as the spectrum peak, if the spectrum peak has a plurality of spectrum peaks, extracting the top m maximum values and the corresponding maximum value points in the E (n) as the spectrum peaks, wherein m is equal to the number of the spectrum peaks.
2. The method for detecting the spectral peak of the underwater acoustic signal according to claim 1, wherein: the moving average filter length in S3 is equal to the short sliding window length in S2.
3. The method for detecting the spectral peak of the underwater acoustic signal according to claim 1, wherein: the length of the gaussian window function in S4 is equal to the length of the short sliding window in S2.
4. The method for detecting the spectral peak of the underwater acoustic signal according to claim 1, wherein: the short sliding window length in S2, the moving average filter length in S3, and the gaussian window function length in S4 are equal.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202111411551.3A CN114152329B (en) | 2021-11-25 | 2021-11-25 | Method for detecting spectral peaks of underwater acoustic signals |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202111411551.3A CN114152329B (en) | 2021-11-25 | 2021-11-25 | Method for detecting spectral peaks of underwater acoustic signals |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN114152329A true CN114152329A (en) | 2022-03-08 |
| CN114152329B CN114152329B (en) | 2023-07-21 |
Family
ID=80457441
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202111411551.3A Active CN114152329B (en) | 2021-11-25 | 2021-11-25 | Method for detecting spectral peaks of underwater acoustic signals |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN114152329B (en) |
Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| FR2618625A1 (en) * | 1987-07-24 | 1989-01-27 | Labo Electronique Physique | DEVICE FOR ENHANCING THE DECODING OF DIGITAL SIGNALS DURING FREQUENCY MODULATION TRANSMISSIONS |
| JP2005099112A (en) * | 2003-09-22 | 2005-04-14 | Nippon Sharyo Seizo Kaisha Ltd | Spectrum peak flattening processing for adaptive control |
| CN101782426A (en) * | 2010-01-29 | 2010-07-21 | 浙江大学 | Detection method of looseness fault vibration of power transformer winding |
| CN101951276A (en) * | 2010-09-30 | 2011-01-19 | 哈尔滨工程大学 | Method for detecting and suppressing Gaussian fitting linear frequency-modulated jamming in direct sequence spread spectrum (DSSS) communication system |
| KR20210016140A (en) * | 2019-07-31 | 2021-02-15 | 삼성전자주식회사 | Spectrum analyzer and method of analyzing spectrum |
| CN112415078A (en) * | 2020-11-18 | 2021-02-26 | 深圳市步锐生物科技有限公司 | Mass spectrum data spectrogram signal calibration method and device |
-
2021
- 2021-11-25 CN CN202111411551.3A patent/CN114152329B/en active Active
Patent Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| FR2618625A1 (en) * | 1987-07-24 | 1989-01-27 | Labo Electronique Physique | DEVICE FOR ENHANCING THE DECODING OF DIGITAL SIGNALS DURING FREQUENCY MODULATION TRANSMISSIONS |
| JP2005099112A (en) * | 2003-09-22 | 2005-04-14 | Nippon Sharyo Seizo Kaisha Ltd | Spectrum peak flattening processing for adaptive control |
| CN101782426A (en) * | 2010-01-29 | 2010-07-21 | 浙江大学 | Detection method of looseness fault vibration of power transformer winding |
| CN101951276A (en) * | 2010-09-30 | 2011-01-19 | 哈尔滨工程大学 | Method for detecting and suppressing Gaussian fitting linear frequency-modulated jamming in direct sequence spread spectrum (DSSS) communication system |
| KR20210016140A (en) * | 2019-07-31 | 2021-02-15 | 삼성전자주식회사 | Spectrum analyzer and method of analyzing spectrum |
| CN112415078A (en) * | 2020-11-18 | 2021-02-26 | 深圳市步锐生物科技有限公司 | Mass spectrum data spectrogram signal calibration method and device |
Also Published As
| Publication number | Publication date |
|---|---|
| CN114152329B (en) | 2023-07-21 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN108831499B (en) | Speech enhancement method using speech existence probability | |
| US11056130B2 (en) | Speech enhancement method and apparatus, device and storage medium | |
| KR100745977B1 (en) | Apparatus and method for voice activity detection | |
| WO2006112831A1 (en) | Method and system for estimating time of arrival of signals using multiple different time scales | |
| CN102546061A (en) | Self-adaptive time-frequency hole detection method based on wavelet transformation | |
| Akay et al. | Broadband interference excision in spread spectrum communication systems via fractional Fourier transform | |
| CN109087657B (en) | Voice enhancement method applied to ultra-short wave radio station | |
| CN103268766A (en) | Method and device for speech enhancement with double microphones | |
| CN108039182B (en) | Voice activation detection method | |
| CN114152329B (en) | Method for detecting spectral peaks of underwater acoustic signals | |
| CN109076038B (en) | Method for estimating parameters of a signal contained in a frequency band | |
| CN107809256A (en) | A kind of shortwave suppressing method under arrowband interference | |
| CN115378776A (en) | MFSK modulation identification method based on cyclic spectrum parameters | |
| CN109660475B (en) | A kind of non-cooperation phase code water sound communication signal autonomous identifying method | |
| US8374229B2 (en) | Method for the detection and generation of a useful signal and associated devices and communications system | |
| CN114650108B (en) | Method and system for detecting signal of transform domain communication system | |
| KR20160116440A (en) | SNR Extimation Apparatus and Method of Voice Recognition System | |
| CN109617839A (en) | A kind of Morse signal detection method based on Kalman filtering algorithm | |
| CN113489552B (en) | Frequency hopping signal detection method based on time-frequency spectrum matrix local variance | |
| WO2002056303A2 (en) | Noise filtering utilizing non-gaussian signal statistics | |
| CN111933169B (en) | Voice noise reduction method for secondarily utilizing voice existence probability | |
| Oliinyk et al. | Center weighted median filter application to time delay estimation in non-Gaussian noise environment | |
| Yeste-Ojeda et al. | Cyclostationarity-based signal separation in interceptors based on a single sensor | |
| CN113593599A (en) | Method for removing noise signal in voice signal | |
| CN109155883A (en) | Noise measuring and noise reduce |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |