CN114842867A - DFT-based audio sinusoidal signal frequency estimation method and system - Google Patents
DFT-based audio sinusoidal signal frequency estimation method and system Download PDFInfo
- Publication number
- CN114842867A CN114842867A CN202210423098.6A CN202210423098A CN114842867A CN 114842867 A CN114842867 A CN 114842867A CN 202210423098 A CN202210423098 A CN 202210423098A CN 114842867 A CN114842867 A CN 114842867A
- Authority
- CN
- China
- Prior art keywords
- spectral line
- frequency
- estimation
- dft
- estimate
- 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.)
- Pending
Links
- 238000000034 method Methods 0.000 title claims abstract description 41
- 230000003595 spectral effect Effects 0.000 claims abstract description 128
- 230000014509 gene expression Effects 0.000 claims description 46
- 238000004364 calculation method Methods 0.000 claims description 26
- 108010076504 Protein Sorting Signals Proteins 0.000 claims description 15
- 238000005070 sampling Methods 0.000 claims description 10
- 238000004422 calculation algorithm Methods 0.000 abstract description 42
- 238000001228 spectrum Methods 0.000 abstract description 12
- 238000010586 diagram Methods 0.000 description 13
- 238000004590 computer program Methods 0.000 description 7
- 238000011161 development Methods 0.000 description 7
- 238000012545 processing Methods 0.000 description 7
- 230000006870 function Effects 0.000 description 4
- 230000005236 sound signal Effects 0.000 description 4
- 230000008569 process Effects 0.000 description 3
- 238000004088 simulation Methods 0.000 description 3
- 230000008859 change Effects 0.000 description 2
- 238000006243 chemical reaction Methods 0.000 description 2
- 238000007476 Maximum Likelihood Methods 0.000 description 1
- 230000004888 barrier function Effects 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000004891 communication Methods 0.000 description 1
- 230000003247 decreasing effect Effects 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 239000000284 extract Substances 0.000 description 1
- 238000004519 manufacturing process Methods 0.000 description 1
- 238000010295 mobile communication Methods 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000003287 optical effect Effects 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 230000009466 transformation Effects 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G10—MUSICAL INSTRUMENTS; ACOUSTICS
- G10L—SPEECH ANALYSIS TECHNIQUES OR SPEECH SYNTHESIS; SPEECH RECOGNITION; SPEECH OR VOICE PROCESSING TECHNIQUES; SPEECH OR AUDIO CODING OR DECODING
- G10L25/00—Speech or voice analysis techniques not restricted to a single one of groups G10L15/00 - G10L21/00
- G10L25/03—Speech or voice analysis techniques not restricted to a single one of groups G10L15/00 - G10L21/00 characterised by the type of extracted parameters
- G10L25/18—Speech or voice analysis techniques not restricted to a single one of groups G10L15/00 - G10L21/00 characterised by the type of extracted parameters the extracted parameters being spectral information of each sub-band
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Spectroscopy & Molecular Physics (AREA)
- Computational Linguistics (AREA)
- Signal Processing (AREA)
- Health & Medical Sciences (AREA)
- Audiology, Speech & Language Pathology (AREA)
- Human Computer Interaction (AREA)
- Acoustics & Sound (AREA)
- Multimedia (AREA)
- Measurement Of Mechanical Vibrations Or Ultrasonic Waves (AREA)
Abstract
The invention discloses a DFT-based audio sinusoidal signal frequency estimation method and system, wherein the method comprises the following steps: firstly, performing spectrum peak search on a DFT magnitude spectrum of a signal to obtain an index value corresponding to a peak spectral line and obtain a frequency rough estimation result; and then, calculating the deviation delta of the frequency rough estimation value through interpolation of left and right spectral lines of the peak spectral line and left and right half spectral lines so as to obtain the complete signal frequency. The method can effectively inhibit white noise, has higher estimation precision and anti-noise performance, and can solve the technical problem of contradiction between algorithm complexity and algorithm estimation precision in the prior art.
Description
Technical Field
The invention belongs to the technical field of signal processing and communication, and particularly relates to a frequency estimation method and system of an audio sinusoidal signal based on DFT (discrete Fourier transform).
Background
The problem of frequency estimation of sinusoidal signals under a noise background is an important research content in the field of digital signal processing, and along with the development of information technology, the frequency estimation method is widely applied to various engineering fields such as mobile communication, radars, audios and videos and the like; an audio signal is a typical signal type suitable for description using a sinusoidal signal model.
At present, many domestic and foreign documents already propose frequency estimation methods for sinusoidal signals, and existing estimation schemes can be divided into two categories, namely a time domain estimation algorithm and a frequency domain estimation algorithm according to the used signal characteristics. The algorithm for estimating the frequency parameter based on the signal time domain characteristic sinusoidal signal mainly comprises a maximum likelihood estimation algorithm, an autocorrelation method, a linear prediction algorithm and the like; the sinusoidal signal estimation algorithm based on the signal frequency domain features usually performs discrete fourier transform on an observed sampled signal, and then extracts parameters such as frequency according to the spectral features of the signal.
Due to the advantages of simple implementation, high calculation efficiency and the like, the DFT estimator becomes a method which is widely concerned and expanded in the field of sinusoidal signal frequency estimation, but inherent defects of spectrum leakage, barrier effect and the like in an incoherent sampling environment can affect the estimation accuracy of the estimator. In the improved frequency estimation algorithm based on DFT, the accuracy of the DFT frequency estimator can be improved to a certain extent by improving measures such as windowing, interpolation and the like; in order to further improve the estimation accuracy, an iterative DFT algorithm is proposed, however, the introduction of iteration generates an additional calculation amount, which causes a huge calculation load and limits the real-time processing capability of the algorithm.
Disclosure of Invention
The present invention is directed to a method and system for frequency estimation of an audio sinusoidal signal based on DFT, so as to solve one or more of the above technical problems. The technical scheme provided by the invention can solve the technical problem of contradiction between algorithm complexity and algorithm estimation precision in the prior art.
In order to achieve the purpose, the invention adopts the following technical scheme:
the invention provides an audio sinusoidal signal frequency estimation method based on DFT, which comprises the following steps:
acquiring a single-tone audio sinusoidal signal to be processed to obtain observation signal sequences of N sample points; wherein the sampling result in the noise-free environment is represented as a signal waveform
A. f and
respectively, the amplitude, frequency and phase of the sinusoidal signal, N-0, 1, …, N-1 being the discrete signal index,
is shown as
v is the digital frequency of the signal, p belongs to {0, 1, …, N-1}, and | delta | less than or equal to 0.5 are respectively the integer part frequency and the decimal part frequency of the digital frequency v;
performing N-point discrete Fourier transform on the observation signal sequences of the N sample points to obtain a transformed signal DFT spectral line sequence X [ k ]]The expression is as follows,
wherein k is 0, 1, …, and N-1 is DFT coefficient index; i is an imaginary symbol;
searching the spectral line with the maximum positioning amplitude of the DFT spectral line sequence X [ k ] as a peak spectral line X [ p ], and taking the index value of the peak spectral line X [ p ] as the estimated value of the integer part frequency p;
calculating said peak spectral line X [ p ]]Is half-indexed interpolated spectral line X [ p + -0.5] in both left and right directions]A spectral line used for determination of estimation is judged based on the spectral line value obtained by calculation; calculating to obtain fractional part frequency based on spectral line used for estimation determined by judgment
An estimated value of (d);
combining the estimated value of the integer part frequency p and the fractional part frequency
Obtaining an estimation result of the digital angular frequency
And finishing the frequency estimation of the audio sinusoidal signal.
In a further development of the method according to the invention, the peak line X [ p ] is calculated]Is half-indexed interpolated spectral line X [ p + -0.5] in both left and right directions]Based on the value ofCalculating the obtained spectral line value, judging and determining the estimated spectral line; calculating to obtain fractional part frequency based on spectral line used for estimation determined by judgment
The step of estimating specifically comprises:
when X [ p-1]]>X[p+0.5]When, X [ p-1]]、X[p-0.5]And X [ p ]]Is the maximum three DFT coefficients, where X [ p-1] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
When X [ p +1]]>X[p-0.5]When, X [ p +1]]、X[p+0.5]And X [ p ]]Is the largest three DFT coefficients, where X [ p +1] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
In other cases, X [ p-0.5]]、X[p+0.5]And X [ p ]]Is the largest three DFT coefficients, where X [ p +0.5] is used]And X [ p-0.5]To obtain an estimate of δ, expressed as
In a further development of the method according to the invention, the peak line X [ p ] is calculated]Is half-indexed interpolated spectral line X [ p + -0.5] in both left and right directions]A spectral line used for determination of estimation is judged based on the spectral line value obtained by calculation; calculating to obtain fractional part frequency based on spectral line used for estimation determined by judgment
The step of estimating specifically comprises:
when X [ p-0.5]>X[p+0.5]When, X [ p-0.5]And X [ p ]]Is the largest two DFT coefficients, where X [ p-0.5] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
When X [ p-0.5]<X[p+0.5]When, X [ p +0.5]And X [ p ]]Is the largest two DFT coefficients, where X [ p +0.5] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
In a further development of the method according to the invention, the peak line X [ p ] is calculated]Is half-indexed interpolated spectral line X [ p + -0.5] in both left and right directions]A spectral line used for determination of estimation is judged based on the spectral line value obtained by calculation; calculating to obtain fractional part frequency based on spectral line used for estimation determined by judgment
The step of estimating specifically comprises:
when X [ p-0.5]>X[p+0.5]When, if
In this case, X [ p-1] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
If it is not
In this case, X [ p-0.5] is used]And X [ p +0.5]To obtain a value of estimated delta, the computational expression being
In other cases, X [ p-0.5] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
When X [ p-0.5]<X[p+0.5]When, if
In this case, X [ p +1] is used]And X [ p ]]To obtain an estimate of deltaThe calculation expression is
If it is not
In this case, X [ p-0.5] is used]And X [ p +0.5]To obtain an estimate of δ, the computational expression is
In other cases, X [ p +0.5] is used in this case]And X [ p ]]To obtain a value of estimated delta, the computational expression being
The invention provides an audio frequency sinusoidal signal frequency estimation system based on DFT, comprising:
the observation signal sequence acquisition module is used for acquiring single-tone audio sinusoidal signals to be processed and acquiring observation signal sequences of N sample points; wherein the sampling result in the noise-free environment is represented as a signal waveform
A. f and
respectively, the amplitude, frequency and phase of the sinusoidal signal, N-0, 1, …, N-1 being the discrete signal index,
is shown as
v is the digital frequency of the signal, p belongs to {0, 1, …, N-1}, and | delta | less than or equal to 0.5 are respectively the integer part frequency and the decimal part frequency of the digital frequency v;
an integer part frequency estimation value acquisition module, configured to perform N-point discrete fourier transform on the observation signal sequence of the N sample points to obtain a transformed signal DFT spectral line sequence X [ k ]]The expression is as follows,
wherein k is 0, 1, …, and N-1 is DFT coefficient index; i is an imaginary symbol; searching the signal DFT spectral line sequence X [ k ]]Locating the maximum spectral line of amplitude as the peak spectral line X [ p ]]The peak spectral line X [ p ]]As an estimate of the integer part frequency p;
a fractional part frequency estimation value acquisition module for calculating the peak spectral line X [ p ]]Is half-indexed interpolated spectral line X [ p + -0.5] in both left and right directions]A spectral line used for determination of estimation is judged based on the spectral line value obtained by calculation; calculating to obtain fractional part frequency based on spectral line used for estimation determined by judgment
An estimated value of (d);
an estimation result acquisition module for combining the estimation value of the integer part frequency p and the fractional part frequency
Obtaining an estimation result of the digital angular frequency
And finishing the frequency estimation of the audio sinusoidal signal.
In a further development of the inventive system, the calculation of the peak spectral line X [ p ] is carried out]Is half-indexed interpolated spectral line X [ p + -0.5] in both left and right directions]A spectral line used for determination of estimation is judged based on the spectral line value obtained by calculation; calculating to obtain fractional part frequency based on spectral line used for estimation determined by judgment
The step of estimating specifically comprises:
when X [ p-1]]>X[p+0.5]When, X [ p-1]]、X[p-0.5]And X [ p ]]Is the maximum three DFT coefficients, where X [ p-1] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
When X [ p +1]]>X[p-0.5]When, X [ p +1]]、X[p+0.5]And X [ p ]]Is the largest three DFT coefficients, where X [ p +1] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
In other cases, X [ p-0.5]]、X[p+0.5]And X [ p ]]Is the maximum three DFT coefficients, where X [ p +0.5] is used]And X [ p-0.5]To obtain an estimate of δ, expressed as
In a further development of the inventive system, the calculation of the peak spectral line X [ p ] is carried out]Is half-indexed interpolated spectral line X [ p + -0.5] in both left and right directions]A spectral line used for determination of estimation is judged based on the spectral line value obtained by calculation; calculating to obtain fractional part frequency based on spectral line used for estimation determined by judgment
The step of estimating specifically comprises:
when X [ p-0.5]>X[p+0.5]When, X [ p-0.5]And X [ p ]]Is the largest two DFT coefficients, where X [ p-0.5] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
When X [ p-0.5]<X[p+0.5]When, X [ p +0.5]And X [ p ]]Is the largest two DFT coefficients, where X [ p +0.5] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
In a further development of the inventive system, the calculation of the peak spectral line X [ p ] is carried out]Is half-indexed interpolated spectral line X [ p + -0.5] in both left and right directions]The value of (a) is,judging and determining a spectral line used for estimation based on the spectral line value obtained by calculation; calculating to obtain fractional part frequency based on spectral line used for estimation determined by judgment
The step of estimating specifically comprises:
when X [ p-0.5]>X[p+0.5]When, if
In this case, X [ p-1] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
If it is not
In this case, X [ p-0.5] is used]And X [ p +0.5]To obtain a value of estimated delta, the computational expression being
In other cases, X [ p-0.5] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
When X [ p-0.5]<X[p+0.5]When, if
In this case, X [ p +1] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
If it is not
In this case, X [ p-0.5] is used]And X [ p +0.5]To obtain an estimate of δ, the computational expression is
OthersIn this case, X [ p +0.5] is used]And X [ p ]]To obtain a value of estimated delta, the computational expression being
Compared with the prior art, the invention has the following beneficial effects:
the invention provides a single-tone sinusoidal signal frequency estimation method based on interpolation DFT, firstly, performing spectrum peak search on a DFT magnitude spectrum of a signal to obtain an index value corresponding to a peak spectral line and obtain a frequency rough estimation result; and then, calculating the deviation delta of the frequency rough estimation value through interpolation of left and right spectral lines and left and right half spectral lines of the peak spectral line so as to obtain complete signal frequency. Specifically, the DFT coefficient with larger amplitude is selected by judgment, so that the result estimated by the estimator is higher in precision when the signal receives noise interference; in addition, the invention uses polynomial fitting, and can realize the utilization of DFT coefficient combination which can not be utilized through a formula; in addition, the invention performs joint complementation on the two algorithms, and can improve the estimation precision on the premise of not increasing the calculation complexity.
Drawings
In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the drawings used in the description of the embodiments or the prior art are briefly introduced below; it is obvious that the drawings in the following description are some embodiments of the invention, and that for a person skilled in the art, other drawings can be derived from them without inventive effort.
FIG. 1 is a flow chart of a DFT-based audio sinusoidal signal frequency estimation method according to an embodiment of the present invention;
FIG. 2 is a schematic diagram of the relationship between DFT coefficients and DTFT transform for monophonic sinusoidal signals;
FIG. 3 is a diagram of the mean square error of the algorithm of the present invention at different decimal frequencies under the condition of 0dB of signal-to-noise ratio;
FIG. 4 is a diagram of the mean square error of the algorithm of the present invention at different decimal frequencies under the condition of a signal-to-noise ratio of 40 dB;
FIG. 5 is a block diagram of the flow of the calculation steps of an algorithm according to an embodiment of the invention;
FIG. 6 is a diagram illustrating the mean square error of the algorithm under different SNR conditions for a monophonic audio signal according to an embodiment of the present invention.
Detailed Description
In order to make the technical solutions of the present invention better understood, the technical solutions in the embodiments of the present invention will be clearly and completely described below with reference to the drawings in the embodiments of the present invention, and it is obvious that the described embodiments are only a part of the embodiments of the present invention, and not all of the embodiments. All other embodiments, which can be derived by a person skilled in the art from the embodiments given herein without making any creative effort, shall fall within the protection scope of the present invention.
It should be noted that the terms "first," "second," and the like in the description and claims of the present invention and in the drawings described above are used for distinguishing between similar elements and not necessarily for describing a particular sequential or chronological order. It is to be understood that the data so used is interchangeable under appropriate circumstances such that the embodiments of the invention described herein are capable of operation in other sequences than those illustrated or described herein. Furthermore, the terms "comprises," "comprising," and "having," and any variations thereof, are intended to cover a non-exclusive inclusion, such that a process, method, system, article, or apparatus that comprises a list of steps or elements is not necessarily limited to those steps or elements expressly listed, but may include other steps or elements not expressly listed or inherent to such process, method, article, or apparatus.
The invention is described in further detail below with reference to the accompanying drawings:
referring to fig. 1, a method for estimating a frequency of a sinusoidal signal based on DFT according to an embodiment of the present invention includes the following steps:
step 1: collecting monophony tones to be processedFrequency sine signal, sampling result in noise-free environment being signal waveform
Wherein A, f and
respectively representing the amplitude, frequency and phase of the sinusoidal signal, N is 0, 1, …, and N-1 is a discrete signal index; for the observed signal sequence of N sample points,
can be expressed as
Wherein v is the digital frequency of the signal, and p is the integer and the decimal part of the digital frequency v respectively, wherein p belongs to {0, 1, …, N-1} and | δ | ≦ 0.5;
step 2: firstly, making N-point discrete Fourier transform on the above-mentioned collected N-point sinusoidal sequence to obtain transformed frequency spectrum Xk,
wherein k is 0, 1, …, and N-1 is DFT coefficient index; i is an imaginary symbol;
and step 3: searching a signal DFT spectral line sequence X [ k ] to locate a spectral line with the maximum amplitude to obtain a rough estimation result of the system frequency, and obtaining an estimation value of an integer part frequency p;
and 4, step 4: calculating the value of an interpolation point X [ p +/-0.5 ] of half index in the left direction and the right direction of the peak spectral line X [ p ];
and 5: determining a point used for estimation based on the determination, and then calculating a fractional part of the estimation
Combining integer part frequenciesThe estimation result of the digital angular frequency can be finally obtained
The step 5 of calculating the estimated value of the fractional part according to the judgment comprises the following steps:
step 5.1: when X [ p-0.5]>X[p+0.5]When, if
In this case, X [ p-1] is used]And X [ p ]]To obtain an estimate of delta, i.e.
Step 5.2: if it is not
In this case, X [ p-0.5] is used]And X [ p +0.5]To obtain a value of estimated delta, i.e.
Step 5.3: at X [ p-0.5]>X[p+0.5]In other cases of time, X [ p-0.5] is used]And X [ p ]]To obtain an estimate of delta, i.e.
Step 5.4: when X [ p-0.5]<X[p+0.5]When, if
In this case, X [ p +1] is used]And X [ p ]]To obtain an estimate of delta, i.e.
Step 5.5: if it is not
In this case, X [ p-0.5] is used]And X [ p +0.5]To obtain an estimate of delta, i.e.
Step 5.6: in other cases, X [ p +0.5] is used in this case]And X [ p ]]To obtain the value of the estimated delta, i.e.
Finally adding the coarse estimation value to obtain the final frequency
The embodiment of the invention discloses a sinusoidal signal frequency estimation method based on DFT, which is suitable for the frequency estimation problem in an audio system. There are now many frequency estimation algorithms based on the DFT domain, which can be mainly classified into two categories, iterative and non-iterative. Most of the existing non-iterative methods have low complexity, but it is difficult to obtain a satisfactory estimation accuracy. The iterative DFT algorithm greatly improves the estimation accuracy, however, the introduction of the iteration generates additional computation. Aiming at the problem, an IpDFT algorithm based on signal adjacent spectral lines is provided. Firstly, performing discrete Fourier transform on a time domain signal, and searching a DFT magnitude spectrum to obtain an index value corresponding to a peak spectral line, namely a rough estimation value p of frequency; and then, calculating the deviation delta of the frequency rough estimation value through interpolation of left and right spectral lines of the peak spectral line and left and right half spectral lines so as to obtain the complete signal frequency. The method provided by the embodiment of the invention can effectively inhibit white noise and has higher estimation precision and anti-noise performance.
Referring to fig. 2 to 6, a method for estimating a frequency of a sinusoidal signal based on DFT according to an embodiment of the present invention includes:
under the noiseless environment, the single tone audio frequency sine signal wave form x [ n ]]Can be expressed in discrete time form as
Wherein A, f and
representing amplitude, frequency and phase, respectively, of a sinusoidal signal, f S Denotes the sampling frequency, N is 0, 1, …, and N-1 is the discrete signal index. For the observed signal sequence of N sample points,
can be expressed as
Where v is the digital frequency of the signal, p ∈ {0, 1, …, N-1} and | δ ≦ 0.5 are the integer and fractional portions of the digital frequency v, respectively.
The following discusses how to estimate the frequency v of a sinusoidal signal by using an interpolation DFT algorithm, including: firstly, N-point discrete Fourier transform is carried out on the acquired N-point sinusoidal sequence to obtain a transformed frequency spectrum X [ k ]. And obtaining an estimated value of the frequency p of the integer part of the signal according to the spectral line sequence after DFT transformation, namely a rough estimation result of the system frequency, which is an index represented by the spectral line with the maximum DFT coefficient. According to the characteristic of DFT, the DFT conversion of complex sinusoidal signals is obtained by translating the result of DTFT conversion of a rectangular window with the length of N on a frequency spectrum and then sampling. FIG. 2 shows the relationship between DTFT and DFT, and it can be seen that the spectrum has a large main lobe and many other side lobes; the amplitude of the main lobe is much larger than that of the side lobe, and for a noise-containing DFT spectrum, better noise interference resistance can be achieved by estimating the frequency by using DFT coefficients with larger amplitude as much as possible, because for the coefficients with larger amplitude, when the noise interference is received, the amplitude of the amplitude change is smaller than the coefficients with smaller amplitude, so that the noise interference resistance is stronger. Therefore, the embodiment of the invention selects the DFT coefficient used for calculating the deviation delta of the frequency rough estimation value by comparing the peak spectral line X [ p ], the left and right spectral lines X [ p-1] and X [ p +1] and the interpolation values X [ p-0.5] and X [ p +0.5] of the left and right half spectral lines, so as to obtain the complete signal frequency.
In the embodiment of the present invention, three points having the largest amplitude among the above five points are obtained. As can be seen from FIG. 2, when X [ p-1]]>X[p+0.5]When, X [ p-1]]，X[p-0.5]，X[p]Is the maximum three DFT coefficients, where X [ p-1] is used]And X [ p ]]To obtain an estimate of δ, expressed as
When X [ p +1]]>X[p-0.5]When, X [ p +1]]，X[p+0.5]，X[p]Is the largest three DFT coefficients, where X [ p +1] is used]And X [ p ]]To obtain an estimate of δ, expressed as
In other cases, X [ p-0.5]]，X[p+0.5]，X[p]Is the largest three DFT coefficients, where X [ p +0.5] is used]And X [ p-0.5]To obtain an estimate of δ, expressed as
Finally adding the rough estimation value to obtain the final frequency
In the above method, the largest DFT coefficient is not used for estimation, because if the largest next largest coefficient is used, the estimated value of δ cannot be obtained from the two coefficients using analytical expressions, where,
however, the ratio of the two coefficients and δ are in one-to-one correspondence, and the functional relationship between the two coefficients is monotonous, so a polynomial fitting method is adopted to describe δ and δ
The relationship between them. The embodiments of the present invention are inventive in that,in order to obtain lower estimator system deviation, the decision is made to fit by using a ninth-order polynomial, and the result of the fitted polynomial is used
To indicate. For delta and
the relationship between the two can be seen by formula
When X [ p-0.5]>X[p+0.5]When, X [ p-0.5]，X[p]Is the largest two DFT coefficients, where X [ p-0.5] is used]And X [ p ]]To obtain an estimate of δ, expressed as
When X [ p-0.5]<X[p+0.5]When, X [ p +0.5]，X[p]Is the largest two DFT coefficients, where X [ p +0.5] is used]And X [ p ]]To obtain an estimate of δ, expressed as
Finally adding the rough estimation value to obtain the final frequency
Illustratively, in order to explore the performance of the algorithm provided by the embodiment of the present invention, the accuracy of the two algorithms is compared after modeling the signal. The signal amplitude A selected by the invention is 1, the phase position is any value, the number of observed sample points is 2048, and the sampling frequency f S At 44100Hz, this signal is a typical high-pitched monophonic audio signal. In order to better evaluate the frequency estimation performance of the method of the invention, a cramer-circle lower bound (CRLB) is introduced as a reference, the CRLB represents the lowest value that the mean square error of the signal parameter estimation can reach, i.e. the optimal performance limit of the estimation algorithm, and for the complex sinusoidal signal model used in the invention,the calculation formula of the frequency estimation CRLB is as follows:
In order to improve the performance better, the invention combines the two algorithms and provides a combined algorithm, and the estimation precision of the algorithm is very high and is superior to that of all the non-iterative DFT domain algorithms at present. The method comprises the following specific steps:
when X [ p-0.5]>X[p+0.5]When, if
In this case, X [ p-1] is used]And X [ p ]]To obtain an estimate of delta, expressed as,
if it is not
In this case, X [ p-0.5] is used]And X [ p +0.5]To obtain a value of the estimated delta, expressed as,
in other cases, X [ p-0.5] is used]And X [ p ]]To obtain an estimate of δ, expressed as
When X [ p-0.5]<X[p+0.5]When, if
In this case, X [ p +1] is used]And X [ p ]]To obtain an estimate of delta, expressed as,
if it is not
In this case, X [ p-0.5] is used]And X [ p +0.5]To obtain an estimate of delta, expressed as,
in other cases, we use X [ p +0.5] at this time]And X [ p ]]To obtain a value of estimated delta, expressed as
Finally adding the rough estimation value to obtain the final frequency
The final algorithm steps are shown in fig. 5.
In combination with the above embodiments, in order to analyze the anti-noise performance of the proposed method more deeply, the monophonic audio signal as mentioned above is used, and the frequency estimation performance of the proposed algorithm is observed along with the change of the SNR of the signal-to-noise ratio, as shown in fig. 6, the mean square error of the frequency estimation of the proposed combined algorithm (combined) of the present invention is approximately linearly decreased along with the increase of the SNR, and the estimation accuracy of the algorithm is improved along with the increase of the SNR; and under the condition of most signal-to-noise ratios, the joint algorithm provided by the invention always has the minimum mean square error and the highest estimation precision in the current non-iterative algorithm. In conclusion, the interpolation DFT-based single-tone sinusoidal signal frequency estimation method provided by the invention has excellent performance, effectively inhibits the influence of white noise, and has the characteristics of simple operation, high estimation precision, strong anti-noise performance and the like.
The following are embodiments of the apparatus of the present invention that may be used to perform embodiments of the method of the present invention. For details of non-careless mistakes in the embodiment of the apparatus, please refer to the embodiment of the method of the present invention.
In another embodiment of the present invention, a DFT-based audio sinusoidal signal frequency estimation system is provided, including:
the observation signal sequence acquisition module is used for acquiring a single-tone audio sinusoidal signal to be processed and acquiring observation signal sequences of N sample points; wherein the sampling result in the noise-free environment is represented as a signal waveform
A. f and
respectively, the amplitude, frequency and phase of the sinusoidal signal, N-0, 1, …, N-1 being the discrete signal index,
is shown as
v is the digital frequency of the signal, p is epsilon {0, 1, …, N-1}, and | delta | is less than or equal to 0.5, and is the integer part frequency and the decimal part frequency of the digital frequency v respectively;
an integer part frequency estimation value acquisition module, configured to perform N-point discrete fourier transform on the observation signal sequence of the N sample points to obtain a transformed signal DFT spectral line sequence X [ k ]]The expression is as follows,
wherein k is 0, 1, …, and N-1 is DFT coefficient index; i is an imaginary symbol; searching the signal DFT spectral line sequence X [ k ]]Locating the maximum spectral line of amplitude as the peak spectral line X [ p ]]The peak value spectrum is obtainedLine X [ p ]]As an estimate of the integer part frequency p;
a fractional part frequency estimation value acquisition module for calculating the peak spectral line X [ p ]]Is half-indexed interpolated spectral line X [ p + -0.5] in both left and right directions]A spectral line used for determination of estimation is judged based on the spectral line value obtained by calculation; calculating to obtain fractional part frequency based on spectral line used for estimation determined by judgment
An estimated value of (d);
an estimation result acquisition module for combining the estimation value of the integer part frequency p and the fractional part frequency
Obtaining an estimation result of the digital angular frequency
And finishing the frequency estimation of the audio sinusoidal signal.
As will be appreciated by one skilled in the art, embodiments of the present application may be provided as a method, system, or computer program product. Accordingly, the present application may take the form of an entirely hardware embodiment, an entirely software embodiment or an embodiment combining software and hardware aspects. Furthermore, the present application may take the form of a computer program product embodied on one or more computer-usable storage media (including, but not limited to, disk storage, CD-ROM, optical storage, and the like) having computer-usable program code embodied therein.
The present application is described with reference to flowchart illustrations and/or block diagrams of methods, apparatus (systems), and computer program products according to embodiments of the application. It will be understood that each flow and/or block of the flow diagrams and/or block diagrams, and combinations of flows and/or blocks in the flow diagrams and/or block diagrams, can be implemented by computer program instructions. These computer program instructions may be provided to a processor of a general purpose computer, special purpose computer, embedded processor, or other programmable data processing apparatus to produce a machine, such that the instructions, which execute via the processor of the computer or other programmable data processing apparatus, create means for implementing the functions specified in the flowchart flow or flows and/or block diagram block or blocks.
These computer program instructions may also be stored in a computer-readable memory that can direct a computer or other programmable data processing apparatus to function in a particular manner, such that the instructions stored in the computer-readable memory produce an article of manufacture including instruction means which implement the function specified in the flowchart flow or flows and/or block diagram block or blocks.
These computer program instructions may also be loaded onto a computer or other programmable data processing apparatus to cause a series of operational steps to be performed on the computer or other programmable apparatus to produce a computer implemented process such that the instructions which execute on the computer or other programmable apparatus provide steps for implementing the functions specified in the flowchart flow or flows and/or block diagram block or blocks.
Finally, it should be noted that: the above embodiments are only for illustrating the technical solutions of the present invention and not for limiting the same, and although the present invention is described in detail with reference to the above embodiments, those of ordinary skill in the art should understand that: modifications and equivalents may be made to the embodiments of the invention without departing from the spirit and scope of the invention, which is to be covered by the claims.
Claims (8)
1. A DFT-based audio sinusoidal signal frequency estimation method is characterized by comprising the following steps:
acquiring a single-tone audio sinusoidal signal to be processed to obtain observation signal sequences of N sample points; wherein the sampling result in the noise-free environment is represented as a signal waveform
A. f and
respectively, the amplitude, frequency and phase of the sinusoidal signal, N-0, 1, …, N-1 being the discrete signal index,
is shown as
v is the digital frequency of the signal, p belongs to {0, 1, …, N-1}, and | delta | less than or equal to 0.5 are respectively the integer part frequency and the decimal part frequency of the digital frequency v;
performing N-point discrete Fourier transform on the observation signal sequences of the N sample points to obtain a transformed signal DFT spectral line sequence X [ k ]]The expression is as follows,
wherein k is 0, 1, …, and N-1 is DFT coefficient index; i is an imaginary symbol;
searching the spectral line with the maximum positioning amplitude of the DFT spectral line sequence X [ k ] as a peak spectral line X [ p ], and taking the index value of the peak spectral line X [ p ] as the estimated value of the integer part frequency p;
calculating said peak spectral line X [ p ]]Is half-indexed interpolated spectral line X [ p + -0.5] in both left and right directions]A spectral line used for determination of estimation is judged based on the spectral line value obtained by calculation; calculating to obtain fractional part frequency based on spectral line used for estimation determined by judgment
An estimated value of (d);
combining the estimated value of the integer part frequency p and the fractional part frequency
Obtaining an estimation result of the digital angular frequency
And finishing the frequency estimation of the audio sinusoidal signal.
2. The DFT-based audio sinusoidal signal frequency estimation method according to claim 1, wherein the calculating the peak spectral line X [ p ]]Is half-indexed interpolated spectral line X [ p + -0.5] in both left and right directions]A spectral line used for determination of estimation is judged based on the spectral line value obtained by calculation; calculating to obtain fractional part frequency based on spectral line used for estimation determined by judgment
The step of estimating specifically comprises:
when X [ p-1]]>X[p+0.5]When, X [ p-1]]、X[p-0.5]And X [ p ]]Is the maximum three DFT coefficients, where X [ p-1] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
When X [ p +1]]>X[p-0.5]When, X [ p +1]]、X[p+0.5]And X [ p ]]Is the largest three DFT coefficients, where X [ p +1] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
In other cases, X [ p-0.5]]、X[p+0.5]And X [ p ]]Is the largest three DFT coefficients, where X [ p +0.5] is used]And X [ p-0.5]To obtain an estimate of δ, expressed as
3. The DFT-based audio sinusoidal signal frequency estimation method according to claim 1, wherein the calculating the peak spectral line X [ p ]]Is half-indexed interpolated spectral line X [ p + -0.5] in both left and right directions]Based on the value ofCalculating the obtained spectral line value, judging and determining the estimated spectral line; calculating to obtain fractional part frequency based on spectral line used for estimation determined by judgment
The step of estimating specifically comprises:
when X [ p-0.5]>X[p+0.5]When, X [ p-0.5]And X [ p ]]Is the largest two DFT coefficients, where X [ p-0.5] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
When X [ p-0.5]<X[p+0.5]When, X [ p +0.5]And X [ p ]]Is the largest two DFT coefficients, where X [ p +0.5] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
4. The DFT-based audio sinusoidal signal frequency estimation method according to claim 1, wherein the calculating the peak spectral line X [ p ]]Is half-indexed interpolated spectral line X [ p + -0.5] in both left and right directions]A spectral line used for determination of estimation is judged based on the spectral line value obtained by calculation; calculating to obtain fractional part frequency based on spectral line used for estimation determined by judgment
The step of estimating specifically comprises:
when X [ p-0.5]>X[p+0.5]When, if
In this case, X [ p-1] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
If it is not
In this case, X [ p-0.5] is used]And X [ p +0.5]To obtain a value of estimated delta, the computational expression being
In other cases, X [ p-0.5] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
When X [ p-0.5]<X[p+0.5]When, if
In this case, X [ p +1] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
If it is not
In this case, X [ p-0.5] is used]And X [ p +0.5]To obtain an estimate of δ, the computational expression is
In other cases, X [ p +0.5] is used in this case]And X [ p ]]To obtain a value of estimated delta, the computational expression being
5. A DFT-based audio sinusoidal signal frequency estimation system, comprising:
the observation signal sequence acquisition module is used for acquiring single-tone audio sinusoidal signals to be processed and acquiring observation signal sequences of N sample points; wherein the sampling result in the noise-free environment is represented as a signal waveform
A. f and
respectively, the amplitude, frequency and phase of the sinusoidal signal, N-0, 1, …, N-1 being the discrete signal index,
is shown as
v is the digital frequency of the signal, p belongs to {0, 1, …, N-1}, and | delta | less than or equal to 0.5 are respectively the integer part frequency and the decimal part frequency of the digital frequency v;
an integer part frequency estimation value acquisition module, configured to perform N-point discrete fourier transform on the observation signal sequence of the N sample points to obtain a transformed signal DFT spectral line sequence X [ k ]]The expression is as follows,
wherein k is 0, 1, …, and N-1 is DFT coefficient index; i is an imaginary symbol; searching the signal DFT spectral line sequence X [ k ]]Locating the maximum spectral line of amplitude as the peak spectral line X [ p ]]The peak spectral line X [ p ]]As an estimate of the integer part frequency p;
a fractional part frequency estimation value acquisition module for calculating the peak spectral line X [ p ]]Is half-indexed interpolated spectral line X [ p + -0.5] in both left and right directions]A spectral line used for determination of estimation is judged based on the spectral line value obtained by calculation; calculating to obtain fractional part frequency based on spectral line used for estimation determined by judgment
An estimated value of (d);
an estimation result acquisition module for combining the estimation value of the integer part frequency p and the fractional part frequency
Obtaining an estimation result of the digital angular frequency
And finishing the frequency estimation of the audio sinusoidal signal.
6. The DFT-based audio sinusoidal signal frequency estimation system according to claim 5, wherein the calculating the peak spectral line X [ p ]]Is half-indexed interpolated spectral line X [ p + -0.5] in both left and right directions]A spectral line used for determination of estimation is judged based on the spectral line value obtained by calculation; calculating to obtain fractional part frequency based on spectral line used for estimation determined by judgment
The step of estimating specifically comprises:
when X [ p-1]]>X[p+0.5]When, X [ p-1]]、X[p-0.5]And X [ p ]]Is the maximum three DFT coefficients, where X [ p-1] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
When X [ p +1]]>X[p-0.5]When, X [ p +1]]、X[p+0.5]And X [ p ]]Is the largest three DFT coefficients, where X [ p +1] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
In other cases, X [ p-0.5]]、X[p+0.5]And X [ p ]]Is the largest three DFT coefficients, where X [ p +0.5] is used]And X [ p-0.5]To obtain an estimate of δ, expressed as
7. According to claim5 the DFT-based frequency estimation system for audio sinusoidal signals, wherein the calculating of the peak spectral line X [ p ]]Is half-indexed interpolated spectral line X [ p + -0.5] in both left and right directions]A spectral line used for determination of estimation is judged based on the spectral line value obtained by calculation; calculating to obtain fractional part frequency based on spectral line used for estimation determined by judgment
The step of estimating specifically comprises:
when X [ p-0.5]>X[p+0.5]When, X [ p-0.5]And X [ p ]]Is the maximum two DFT coefficients, where X [ p-0.5] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
When X [ p-0.5]<X[p+0.5]When, X [ p +0.5]And X [ p ]]Is the maximum two DFT coefficients, where X [ p +0.5] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
8. The DFT-based audio sinusoidal signal frequency estimation system according to claim 5, wherein the calculating the peak spectral line X [ p ]]Is half-indexed interpolated spectral line X [ p + -0.5] in both left and right directions]A spectral line used for determination of estimation is judged based on the spectral line value obtained by calculation; calculating to obtain fractional part frequency based on the estimated spectral line determined by judgment
The step of estimating specifically comprises:
when X [ p-0.5]>X[p+0.5]When, if
In this case, X [ p-1] is used]And X [ p ]]To obtain an estimate of deltaThe calculation expression is
If it is used
In this case, X [ p-0.5] is used]And X [ p +0.5]To obtain a value of estimated delta, the computational expression being
In other cases, X [ p-0.5] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
When X [ p-0.5]<X[p+0.5]When, if
In this case, X [ p +1] is used]And X [ p ]]To obtain an estimate of δ, the computational expression is
If it is not
In this case, X [ p-0.5] is used]And X [ p +0.5]To obtain an estimate of δ, the computational expression is
In other cases, X [ p +0.5] is used in this case]And X [ p ]]To obtain a value of estimated delta, the computational expression being
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202210423098.6A CN114842867A (en) | 2022-04-21 | 2022-04-21 | DFT-based audio sinusoidal signal frequency estimation method and system |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202210423098.6A CN114842867A (en) | 2022-04-21 | 2022-04-21 | DFT-based audio sinusoidal signal frequency estimation method and system |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN114842867A true CN114842867A (en) | 2022-08-02 |
Family
ID=82566206
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202210423098.6A Pending CN114842867A (en) | 2022-04-21 | 2022-04-21 | DFT-based audio sinusoidal signal frequency estimation method and system |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN114842867A (en) |
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN116125138A (en) * | 2023-04-17 | 2023-05-16 | 湖南工商大学 | Method and device for rapidly estimating frequency of sinusoidal signal based on rotation adjustment |
| CN116257730A (en) * | 2023-05-08 | 2023-06-13 | 成都戎星科技有限公司 | Method for realizing frequency offset tracking based on FPGA |
-
2022
- 2022-04-21 CN CN202210423098.6A patent/CN114842867A/en active Pending
Cited By (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN116125138A (en) * | 2023-04-17 | 2023-05-16 | 湖南工商大学 | Method and device for rapidly estimating frequency of sinusoidal signal based on rotation adjustment |
| CN116257730A (en) * | 2023-05-08 | 2023-06-13 | 成都戎星科技有限公司 | Method for realizing frequency offset tracking based on FPGA |
| CN116257730B (en) * | 2023-05-08 | 2023-08-01 | 成都戎星科技有限公司 | Method for realizing frequency offset tracking based on FPGA |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN110333389B (en) | Sinusoidal signal frequency estimation method based on interpolation DFT | |
| CN108037361B (en) | High-precision harmonic parameter estimation method based on sliding window DFT | |
| CN114842867A (en) | DFT-based audio sinusoidal signal frequency estimation method and system | |
| CN106597408B (en) | High-order PPS signal parameter estimation method based on time-frequency analysis and instantaneous frequency curve fitting | |
| CN112395992A (en) | Electric power harmonic signal denoising method based on improved wavelet threshold | |
| CN112446323A (en) | HHT harmonic analysis method based on improved EMD modal aliasing and endpoint effect | |
| Huibin et al. | Energy based signal parameter estimation method and a comparative study of different frequency estimators | |
| Cheng et al. | A new time delay estimator based on ETDE | |
| Belega et al. | Amplitude and phase estimation of real-valued sine wave via frequency-domain linear least-squares algorithms | |
| Tamim et al. | Techniques for optimization in time delay estimation from cross correlation function | |
| CN116502042A (en) | Power quality disturbance denoising method based on variational modal decomposition and improved wavelet threshold | |
| CN111046791A (en) | Current signal filtering and denoising method based on generalized S transform containing variable factors | |
| Quinn | Recent advances in rapid frequency estimation | |
| CN114785379A (en) | Underwater sound JANUS signal parameter estimation method and system | |
| CN108933746B (en) | Multi-tone interference parameter estimation method based on three-level iteration | |
| Djurovic | Estimation of the sinusoidal signal frequency based on the marginal median DFT | |
| KR102304464B1 (en) | Method And Apparatus for Providing Interpolation for Frequency Offset Estimation | |
| CN112035790A (en) | Method for estimating frequency of positioning signal between wells | |
| CN110112757B (en) | Low-frequency oscillation analysis method based on SURE wavelet denoising and improved HHT | |
| CN112883787B (en) | Short sample low-frequency sinusoidal signal parameter estimation method based on spectrum matching | |
| CN111917676B (en) | Linear frequency modulation interference cancellation method | |
| CN109117698B (en) | Noise background estimation method based on minimum mean square error criterion | |
| CN113468474A (en) | Power grid frequency estimation method based on root Mini-Norm | |
| Zhang et al. | A hybrid phase-based single frequency estimator | |
| CN115310478A (en) | Method for extracting effective Raman spectrum of natural gas based on wavelet optimal decomposition layer number |
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 |