CN112964931B - Non-ideal multi-damping harmonic signal parameter measurement method based on two-channel undersampling - Google Patents

Non-ideal multi-damping harmonic signal parameter measurement method based on two-channel undersampling Download PDF

Info

Publication number
CN112964931B
CN112964931B CN202110103500.8A CN202110103500A CN112964931B CN 112964931 B CN112964931 B CN 112964931B CN 202110103500 A CN202110103500 A CN 202110103500A CN 112964931 B CN112964931 B CN 112964931B
Authority
CN
China
Prior art keywords
signal
sampling
frequency
channel
parameter
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.)
Active
Application number
CN202110103500.8A
Other languages
Chinese (zh)
Other versions
CN112964931A (en
Inventor
黄国兴
倪安
卢为党
彭宏
张昱
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Zhejiang University of Technology ZJUT
Original Assignee
Zhejiang University of Technology ZJUT
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Zhejiang University of Technology ZJUT filed Critical Zhejiang University of Technology ZJUT
Priority to CN202110103500.8A priority Critical patent/CN112964931B/en
Publication of CN112964931A publication Critical patent/CN112964931A/en
Application granted granted Critical
Publication of CN112964931B publication Critical patent/CN112964931B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • 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

Landscapes

  • Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • General Physics & Mathematics (AREA)
  • Magnetic Resonance Imaging Apparatus (AREA)

Abstract

A non-ideal multi-damping harmonic signal parameter measuring method based on double-channel undersampling comprises two sampling channels, in order to prevent frequency mixing and image frequency blurring, the two sampling channels are connected through a feedback channel, signal parameters can be recovered from non-ideal K harmonic MEDS signals through only 4K samples, and the method is also suitable for K unknown conditions. Compared with the method in the prior art, the method has higher reconstruction precision and robustness to noise.

Description

Non-ideal multi-damping harmonic signal parameter measurement method based on double-channel undersampling
Technical Field
The invention belongs to the technical field of signal processing, and particularly relates to a non-ideal multi-damping harmonic signal parameter measurement method based on dual-channel undersampling.
Background
Parameter estimation of multi-Damped harmonic signals (MEDS) is a problem that often occurs in the fields of magnetic resonance imaging, linear system recognition, speech analysis, and transient analysis. A related problem includes pole-zero modeling of empirically generated time series data, since the impulse response of a linear stable rational model is composed of exponentially damped sinusoids. A number of researchers have proposed methods for parameter estimation based on nyquist sampling theory. Fotinea et al propose a state space approach to spectral estimation that performs factor 2 decimation while utilizing all available data sets. Chen et al propose a subspace-based parameter estimator that is efficient and accurate in estimating the parameters of an exponentially decaying sinusoidal signal sum. However, according to the nyquist sampling theorem, the above system needs to satisfy the requirement that the sampling rate is twice the signal bandwidth. Therefore, for wideband signals, these systems face a common problem in that not only a high sampling rate but also a complicated process is required to complete parameter estimation.
The continuous-time med signal is generally composed of a linear combination of a plurality of damped complex exponential signals, but since the actual signal cannot be perfectly matched with the signal, the present invention takes into account the problem of model matching errors and reduces the influence of the errors by an optimization algorithm. Also, due to the periodicity of the trigonometric function, the under-nyquist sampling results in frequency aliasing and image frequency aliasing.
In recent years, a parameter measurement method based on sub-nyquist sampling has received increasing attention. However, when measuring parameters using the under-nyquist sampling method, a key problem to be solved is frequency aliasing. Scholars of r.venkataramani and a.bourdoux, et al, propose multi-rate asynchronous undersnqquist sampling schemes to solve the frequency aliasing problem, which methods typically limit the number of frequency components, also depending on the number of channels. However, these methods often require random sampling, resulting in a complex hardware structure, which greatly reduces the practicality. Furthermore, the estimation accuracy depends on the spectral grid density, which is generally not high in view of computational complexity. In the systems of p.pal and s.qin et al, a two-channel sampling scheme of a relatively prime undersampling rate is employed to estimate the line spectrum. Huang et al, however, show that frequency estimates sometimes cannot be uniquely determined. To solve this problem, s.huang et al propose a three-channel under-Nyquist scheme with pairwise co-prime undersampling rate to ensure that the frequency aliasing problem can be successfully solved. However, these methods require a large number of samples to determine the correct frequency, and these systems require that the number of frequencies be known, and that they are not specifically designed for med signals. Therefore, the measurement problem of the MEDS signal parameters needs to be solved urgently.
Disclosure of Invention
In order to overcome the defects of the prior art, the invention provides a non-ideal multi-damping harmonic signal parameter measurement method based on dual-channel undersampling, which aims at the problems of frequency aliasing and image frequency aliasing and parameter measurement caused by the undersnyquist sampling of a non-ideal MEDS signal; then, an intelligent optimization method is also provided to solve the problems of signal non-ideality and insufficient prior information.
The technical scheme adopted by the invention for solving the technical problems is as follows:
a non-ideal multi-damping harmonic signal parameter measurement method based on dual-channel undersampling comprises the following steps:
firstly, initializing to generate a non-ideal MEDS signal, wherein the signal has the following form:
Figure BDA0002916479670000031
wherein T ∈ (0, T), ck(K-1, 2, … K) is the complex amplitude, sk＝rk+j2πfkIs a complex frequency, K is 1,2, … K, rkIs the damping factor, fkThe method comprises the following steps that pole frequency is adopted, K is the number of unknown frequency components, sigma (t) is a model matching error signal, and the model matching error signal sigma (t) is introduced into a model because a real signal and an MEDS signal cannot be completely matched;
step two, controlling a signal control switch so as to shunt the MEDS signals, firstly, dialing the control switch to a branch where a main sampling channel is located through the control signal, sampling the obtained shunt MEDS signals by the main sampling channel at the moment, then dialing the control switch to a branch where an auxiliary sampling channel is located through the control signal, and sampling the obtained shunt MEDS signals by the auxiliary sampling channel at the moment;
step three, the control signal dials the switch to the branch where the main sampling channel is located, the main sampling channel samples the obtained shunt MEDS signal, the sampling clock CLK1 uniformly samples the shunt signal, and the sampling rate is fsThe sampled values are represented as follows:
Figure BDA0002916479670000032
wherein K is the unknown number of components, ckWhich is indicative of the complex amplitude of the signal,
Figure BDA0002916479670000033
sk＝rk+j2πfkdenotes the complex frequency, rkDenotes the damping factor, fkRepresenting the pole frequency, m ∈ Z+，σ[n]Represents a sample value after σ (t) is sampled by the clock CLK 1;
and step four, establishing an optimization model to reduce the influence of model matching errors. Since the number of components K is unknown and the model matching error sigma (t) ≠ 0, the optimal component K and the optimal parameter are found by minimizing the energy of sigma (t)
Figure BDA0002916479670000034
Step five, sampling values x [ n ]]Input to the optimal estimation algorithm of the number of amplitude and frequency components to obtain the optimal K and
Figure BDA0002916479670000035
p (P ≦ K) indicates a valid estimate of the threshold exceeded
Figure BDA0002916479670000041
The number of (2);
step six, estimating the parameter
Figure BDA0002916479670000042
Input to a feedback sample rate generator module to calculate the sample rate f 'of the secondary sample channel'sAnd f'sSatisfies the following formula:
Figure BDA0002916479670000043
wherein, the first and the second end of the pipe are connected with each other,
Figure BDA0002916479670000044
a,b∈{0,1,…,K-1}，m∈{-(2Q-2),-(2Q-3),…,2Q-2}，
Figure BDA0002916479670000045
fsis the sampling rate, 0 < fs＜fmax，fmaxFor signal frequency upper limit, when the sampling value is interfered by noise, in order to improve the robustness of the method, the sampling rate f 'generated by the sampling rate generator module is fed back'sThe following equation needs to be satisfied:
Figure BDA0002916479670000046
wherein rem { (. cndot.), 1} represents the remainder of (·) divided by 1, and ε represents a threshold determined by the noise strength;
and step seven, the control signal dials the switch to the auxiliary sampling channel, and the auxiliary sampling channel samples the obtained shunt signal. The signal is uniformly sampled by a sampling clock CLK2 at a sampling rate of f'sThe sampled values are represented as follows:
Figure BDA0002916479670000047
wherein, c'kIs a function of the complex amplitude of the signal,
Figure BDA0002916479670000048
T′s＝1/f′sis the sampling period, σ]n′]＝σ[n′T′s]For the sample values of the model matching error signal, N 'is 0,1, …, N' -1, and then an optimal estimation algorithm for the number of amplitude and frequency components is run, the sample values x 'N']By the algorithm, an estimated value is obtained
Figure BDA0002916479670000049
Step eight, estimating the parameter
Figure BDA00029164796700000410
And
Figure BDA00029164796700000411
input into a joint estimation parameter algorithm for complex frequency and amplitude
Figure BDA00029164796700000412
And
Figure BDA00029164796700000413
and (4) obtaining a joint estimation.
Further, the process of the fourth step is as follows:
step 4.1, the established optimization model is expressed as follows:
Figure BDA0002916479670000051
due to the fact that
Figure BDA0002916479670000052
(3) The conversion is:
Figure BDA0002916479670000053
step 4.2, because the input non-ideal med signal x (t) is unknown, that is, the formula (4) cannot be solved, and x (t) is replaced by the sampling value obtained in step (2), then (4) is converted into the following formula:
Figure BDA0002916479670000054
wherein, the first and the second end of the pipe are connected with each other,
Figure BDA0002916479670000055
number of frequency components K and parameter
Figure BDA0002916479670000056
Are unknown, then (5) is described by the following formula:
Figure BDA0002916479670000057
wherein, the first and the second end of the pipe are connected with each other,
Figure BDA0002916479670000058
(7) referred to as the optimization objective function.
Still further, the process of the fifth step is as follows:
step 5.1, initialization, first, obtaining the sampling value x [ n ] from the parameter measurement system](N is 0,1, …, N-1), setting the number of frequency components K' to 0, setting the optimum number of frequency components K to 0, and setting the estimation parameter PbestSetting the optimum value of equation (7) to fbest；
Step 5.2, updating the number K' +1 of the frequency components;
step 5.3, using the sampling value x [ n ]]And solving (2) by a spectrum estimation method to obtain the estimated parameters
Figure BDA0002916479670000059
As an initial bit for the next step;
step 5.4, the optimization problem of (6) is solved by the optimization algorithm, and P is setoptIs a global optimum position;
step 5.5, updating the optimal solution
Figure BDA00029164796700000510
Step 5.6, the result is output, if
Figure BDA0002916479670000061
And (4) outputting the solution to be the optimal solution, otherwise, returning to the step 5.3 and repeating the steps.
The process of the step eight is as follows:
step 8.1, estimated normalized frequency due to the periodicity of the trigonometric function
Figure BDA0002916479670000062
And
Figure BDA0002916479670000063
by a parameter of 2 pi mkNamely:
Figure BDA0002916479670000064
further, the following formula is obtained:
Figure BDA0002916479670000065
wherein the angle (·) is the argument of complex (·), and is not less than 0 ≦ angle (·) < 2 π, mkE.g., Z, K e.g., {1,2, …, K }. Because of
Figure BDA0002916479670000066
fmaxFor the upper frequency limit of the MEDS signal, the following equation is obtained:
Figure BDA0002916479670000067
note book
Figure BDA0002916479670000068
Then (13) is expressed as:
0≤mk≤(Q-1) (14)
step 8.2, according to (12), from
Figure BDA0002916479670000069
The argument principal value of (a) is obtained as a frequency parameter solution set F and expressed as follows:
Figure BDA00029164796700000610
wherein, fsIs the sampling rate of the main sampling channel, the angle is the amplitude main value of the complex number, the angle is more than or equal to 0 and less than 2 pi, mp∈Z，-(Q-1)≤nk≤Q-1，nkE.g. Z, thus, the true frequency
Figure BDA00029164796700000611
While obtaining a damping factorThe estimated values are:
Figure BDA00029164796700000612
step 8.3, in the subsampling channel, according to (12), from
Figure BDA00029164796700000613
The argument of (c) is obtained as a frequency parameter solution set F' and expressed in the form:
Figure BDA0002916479670000071
wherein, f'sIs the sampling rate of the secondary sampling channel, the angle is the amplitude primary value of a complex number, and is more than or equal to 0 and less than 2 pi, m'p′∈Z，-(Q′-1)≤n′k≤Q′-1，n′kE.g. Z, thus, the true frequency
Figure BDA0002916479670000072
Meanwhile, the estimated value of the damping factor is obtained as follows:
Figure BDA0002916479670000073
thus, the frequency parameter is obtained as shown in the following formula:
Figure BDA0002916479670000074
and the damping factor parameter is shown as follows:
Figure BDA0002916479670000075
thereby obtaining a complex frequency skAs shown in the following formula:
Figure BDA0002916479670000076
step 8.4, once the frequency is found
Figure BDA0002916479670000077
We can calculate the parameters
Figure BDA0002916479670000078
And
Figure BDA0002916479670000079
first, by selecting
Figure BDA00029164796700000710
And
Figure BDA00029164796700000711
to determine parameters
Figure BDA00029164796700000712
Assume other elements as
Figure BDA00029164796700000713
And
Figure BDA00029164796700000714
then, the parameters
Figure BDA00029164796700000715
And
Figure BDA00029164796700000716
from the set respectively
Figure BDA00029164796700000717
And
Figure BDA00029164796700000718
is obtained, therefore, parameter ckEstimated as:
Figure BDA00029164796700000719
the invention has the following beneficial effects: obtaining the parameter complex frequency s of the MEDS signal to be measured by using the optimal estimation algorithm and the joint estimation parameter algorithm of the number of the amplitude and frequency componentskAnd complex amplitude ckA value of (d); the method has high reconstruction precision and robustness to noise.
Drawings
FIG. 1 is a system structure diagram of a non-ideal multi-damping harmonic signal parameter measurement method based on two-channel undersampling.
Fig. 2 is a flow chart of an optimal estimation algorithm for the number of amplitude and frequency components.
Fig. 3 is a result of parameter recovery when the frequency is not aliased.
Fig. 4 is a result of parameter recovery at the time of frequency aliasing.
Figure 5 is a comparison of parameter recovery performance as signal-to-noise ratio increases.
Detailed Description
The invention is further described below with reference to the accompanying drawings.
Referring to fig. 1 to 5, a method for measuring non-ideal multi-damping harmonic signal parameters based on dual-channel undersampling comprises the following steps:
step one, initializing to generate a non-ideal MEDS signal, wherein the signal has the following form:
Figure BDA0002916479670000081
wherein T ∈ (0, T), ck(K-1, 2, … K) is the complex amplitude, sk＝rk+j2πfkIs a complex frequency, K is 1,2, … K, rkIs the damping factor, fkIs the pole frequency, K is the unknown number of frequency components, and σ (t) is the model match error signal. Since the true signal and the med signal cannot be completely matched, we introduce a model matching error signal σ (t) in the model;
step two, controlling a signal control switch so as to shunt the MEDS signals, firstly, dialing the control switch to a branch where a main sampling channel is located through the control signal, at the moment, sampling the obtained shunt MEDS signals through the main sampling channel, then dialing the control switch to a branch where an auxiliary sampling channel is located through the control signal, and at the moment, sampling the obtained shunt MEDS signals through the auxiliary sampling channel;
step three, the control signal dials the switch to the branch where the main sampling channel is located, the main sampling channel samples the obtained shunt MEDS signal, the sample clock CLK1 uniformly samples the shunt signal, and the sampling rate is fsThe sample values are represented as follows:
Figure BDA0002916479670000091
wherein K is the unknown number of components, ckWhich is indicative of the complex amplitude of the signal,
Figure BDA0002916479670000092
sk＝rk+j2πfkrepresenting complex frequency, rkDenotes the damping factor, fkRepresenting the pole frequency, m ∈ Z+，σ[n]Represents a sample value after σ (t) is sampled by the clock CLK 1;
step four, establishing an optimization model, reducing the influence of model matching errors, and finding the optimal component K and the optimal parameter by minimizing the energy of sigma (t) because the component number K is unknown and the model matching error sigma (t) is not equal to 0
Figure BDA0002916479670000093
The process is as follows:
step 4.1, the established optimization model is expressed as follows:
Figure BDA0002916479670000094
due to the fact that
Figure BDA0002916479670000095
(3) To convert to:
Figure BDA0002916479670000096
step 4.2, because the input non-ideal med signal x (t) is unknown, that is, the formula (4) cannot be solved, and x (t) is replaced by the sampling value obtained in step (2), then (4) is converted into the following formula:
Figure BDA0002916479670000097
wherein the content of the first and second substances,
Figure BDA0002916479670000098
number of frequency components K and parameter
Figure BDA0002916479670000099
Are unknown, then (5) is described by the following formula:
Figure BDA00029164796700000910
wherein the content of the first and second substances,
Figure BDA00029164796700000911
(7) referred to as an optimization objective function;
step five, sampling value x [ n ]]Input to the optimal estimation algorithm of the number of amplitude and frequency components, thereby obtaining the optimal K and
Figure BDA0002916479670000101
p (P ≦ K) represents a valid estimate of the threshold exceeded
Figure BDA0002916479670000102
The process is as follows:
step 5.1, initialization, first, the sampling values can be obtained from the parameter measuring system shown in FIG. 1x[n](N-0, 1, …, N-1), setting the number of frequency components K' to 0, setting the optimum number of frequency components K to 0, and setting the estimation parameter PbestSetting the optimum value of equation (7) to fbest；
Step 5.2, updating the number K '═ K' +1 of the frequency components;
step 5.3, using the sampling value x [ n ]]And solving (2) by a spectrum estimation method to obtain the estimation parameters
Figure BDA0002916479670000103
As an initial bit for the next step;
step 5.4, the optimization problem of (6) is solved by the provided optimization algorithm, and the P is assumedoptIs a global optimum position;
step 5.5, update the optimal solution
Figure BDA0002916479670000104
Step 5.6, the result is output, if
Figure BDA0002916479670000105
The output is the optimal solution, otherwise, the step 5.3 is returned to, and the steps are repeated;
step six, estimating the parameters
Figure BDA0002916479670000106
Input to a feedback sample rate generator module that calculates the sample rate f 'of the secondary sampling channel'sAnd f'sSatisfies the following formula:
Figure BDA0002916479670000107
wherein, the first and the second end of the pipe are connected with each other,
Figure BDA0002916479670000108
a,b∈{0,1,…,K-1}，m∈{-(2Q-2),-(2Q-3),…,2Q-2}，
Figure BDA0002916479670000109
fsis the sampling rate, 0 < fs＜fmax，fmaxFor the upper limit of the signal frequency, when the sampling value is interfered by noise, in order to improve the robustness of the method, the sampling rate f 'generated by the sampling rate generator module is fed back'sThe following equation needs to be satisfied:
Figure BDA00029164796700001010
wherein rem { (. 1) } represents a remainder of division of (-) by 1, and ε represents a threshold determined by noise intensity;
step seven, the switch is set by the control signal to be in a sub-sampling channel, the sub-sampling channel samples the obtained shunt signal, a sampling clock CLK2 samples the shunt signal uniformly, and the sampling rate is f'sThe sampled values are represented as follows:
Figure BDA0002916479670000111
wherein, c'kIs a function of the complex amplitude of the signal,
Figure BDA0002916479670000112
T′s＝1/f′sis the sampling period, σ [ n']＝σ[n′T′s]For model matching the sampled values of the error signal, N 'is 0,1, …, N' -1, and then an optimal estimation algorithm of the number of amplitude and frequency components is run, the sampled values x 'N']Through the algorithm, an estimated value can be obtained
Figure BDA0002916479670000113
Step eight, estimating the parameter
Figure BDA0002916479670000114
And
Figure BDA0002916479670000115
is input into a joint estimation parameter algorithm,complex frequency and complex amplitude can be used
Figure BDA0002916479670000116
And
Figure BDA0002916479670000117
and (4) obtaining the joint estimation. The algorithm process is as follows:
step 8.1, estimated normalized frequency due to the periodicity of the trigonometric function
Figure BDA0002916479670000118
And with
Figure BDA0002916479670000119
By 2 pi mkNamely:
Figure BDA00029164796700001110
further, the following formula is obtained:
Figure BDA00029164796700001111
wherein, the angle is the main value of the complex angle, and is more than or equal to 0 and less than 2 pi, mkE Z, K e {1,2, …, K }. Because of
Figure BDA00029164796700001112
fmaxFor the upper frequency limit of the MEDS signal, the following equation is obtained:
Figure BDA00029164796700001113
note book
Figure BDA0002916479670000121
Then (13) is expressed as:
0≤mk≤(Q-1) (14)
step 8.2, rootingAccording to (12), from
Figure BDA0002916479670000122
The argument principal value of (a) is obtained as a frequency parameter solution set F and expressed as follows:
Figure BDA0002916479670000123
wherein, fsIs the sampling rate of a main sampling channel, and the angle is the amplitude main value of a complex number, is more than or equal to 0 and less than 2 pi, mp∈Z，-(Q-1)≤nk≤Q-1，nkE.g. Z, thus, the true frequency
Figure BDA0002916479670000124
Meanwhile, the estimated value of the obtained damping factor is as follows:
Figure BDA0002916479670000125
step 8.3, in the sub-sampling channel, according to (12), from
Figure BDA0002916479670000126
The argument of (c) is obtained as a frequency parameter solution set F' and expressed in the form:
Figure BDA0002916479670000127
wherein, f'sIs the sampling rate of the secondary sampling channel, the angle is the amplitude primary value of a complex number, and is more than or equal to 0 and less than 2 pi, m'p′∈Z，-(Q′-1)≤n′k≤Q′-1，n′kE.g. Z, so the true frequency
Figure BDA0002916479670000128
Meanwhile, the estimated value of the obtained damping factor is as follows:
Figure BDA0002916479670000129
thus, the frequency parameter is obtained as shown in the following formula:
Figure BDA00029164796700001210
and the damping factor parameter is as follows:
Figure BDA00029164796700001211
thereby obtaining a complex frequency skAs shown in the following formula:
Figure BDA0002916479670000131
step 8.4, once the frequency is found
Figure BDA0002916479670000132
Calculating parameters
Figure BDA0002916479670000133
And
Figure BDA0002916479670000134
first, by selecting
Figure BDA0002916479670000135
And
Figure BDA0002916479670000136
to determine parameters
Figure BDA0002916479670000137
Let other elements be
Figure BDA0002916479670000138
And
Figure BDA0002916479670000139
then, the parameters
Figure BDA00029164796700001310
And
Figure BDA00029164796700001311
respectively from the set
Figure BDA00029164796700001312
And
Figure BDA00029164796700001313
is obtained, therefore, parameter ckEstimated as:
Figure BDA00029164796700001314
example (c): in a first experiment, it was verified that the method of the present invention can reconstruct the original signal with a small number of sample values under a noise-free condition, and compare it with a clock-interleaved sampling system. In this experiment, the number of frequency components K is set to 5, and the frequency f is set tok＝[0.4,3,5.2,6.8,9.1]KHz. The method of the invention can recover the parameters of the original signal by using 4K to 20 sampling values, and the clock staggered sampling system can also accurately recover the parameters of the original signal, and the recovery result is shown in figure 3.
To further verify the validity of the method of the invention, we set the frequency fk＝[0.4,3,5.5,6.8,9.1]Due to the fact that
Figure BDA00029164796700001315
This will cause aliasing at the frequencies of 3KHz and 5.5KHz, and it can be seen from fig. 4 that the feedback sampling system provided by the present invention can accurately recover the parameters of the original signal, whereas the clock interleaving sampling system cannot accurately recover the parameters of the original signal.
In a second experiment, we compared the method of the invention with several other methods in the presence of noise, includingA clock interleaved sampling system, a time delay sampling system and a three-channel co-prime sampling system. In the experiment process, we set the frequency component number K to 5, and the frequency fkRandomly distributing the signals on (0,10) KHz, and then adding Gaussian noise into the MEDS signals to be measured, wherein the signal-to-noise ratio SNR of the noise is defined as follows:
Figure BDA00029164796700001316
wherein, PsignalAnd PnoiseRespectively representing the power of the signal and the noise, and the SNR ranges from [ -20,50 [)]dB, the number of samples Num obtained by all systems is 300. The experiments are carried out for 2000 times totally, the experimental result is shown in figure 5, and it can be seen from figure 5 that the performance difference between the clock staggered sampling system and the time delay sampling system is not large, when the SNR is less than or equal to 25dB, the performance of the method is similar to that of a three-channel co-prime sampling system, but the method has higher accuracy than other systems along with the increase of the SNR. Experimental results show that the method has better robustness.
The embodiments described in this specification are merely illustrative of implementations of the inventive concepts, which are intended for purposes of illustration only. The scope of the present invention should not be construed as being limited to the particular forms set forth in the embodiments, but is to be accorded the widest scope consistent with the principles and equivalents thereof as contemplated by those skilled in the art.

Claims (4)

1. A non-ideal multi-damping harmonic signal parameter measurement method based on double-channel undersampling is characterized by comprising the following steps:
firstly, initializing to generate a non-ideal MEDS signal, wherein the signal has the following form:
Figure FDA0003612457410000011
wherein T is belonged to (0, T), ckIs the complex amplitude, sk＝rk+j2πfkIs the complex frequency, K is 1,2, … K, rkIs a damping factor, fkThe method comprises the following steps that pole frequency is adopted, K is the number of unknown frequency components, sigma (t) is a model matching error signal, and the model matching error signal sigma (t) is introduced into a model because a real signal and an MEDS signal cannot be completely matched;
step two, controlling a signal control switch so as to shunt the MEDS signals, firstly, dialing the control switch to a branch where a main sampling channel is located through the control signal, sampling the obtained shunt MEDS signals by the main sampling channel at the moment, then dialing the control switch to a branch where an auxiliary sampling channel is located through the control signal, and sampling the obtained shunt MEDS signals by the auxiliary sampling channel at the moment;
step three, the control signal dials the switch to the branch where the main sampling channel is located, the main sampling channel samples the obtained shunt MEDS signal, the sampling clock CLK1 uniformly samples the shunt signal, and the sampling rate is fsThe sample values are represented as follows:
Figure FDA0003612457410000012
where K is the unknown number of components, ckWhich is indicative of the complex amplitude of the signal,
Figure FDA0003612457410000013
sk＝rk+j2πfkrepresenting complex frequency, rkRepresenting the damping factor, fkRepresenting the pole frequency, n ∈ Z+，σ[n]Represents the sample value after σ (t) is sampled by the clock CLK 1;
step four, establishing an optimization model, reducing the influence of model matching errors, and finding the optimal component K and the optimal parameter by minimizing the energy of sigma (t) because the component number K is unknown and the model matching error sigma (t) is not equal to 0
Figure FDA0003612457410000014
Step five, mixingSampling value x [ n ]]Input to the optimal estimation algorithm of the number of amplitude and frequency components to obtain the optimal K and
Figure FDA0003612457410000015
p represents a valid estimate of the exceeding of a threshold
Figure FDA0003612457410000016
P is less than or equal to K;
step six, estimating the parameter
Figure FDA0003612457410000017
Input to a feedback sample rate generator module calculating the sample rate f 'of the secondary sample channel'sAnd f'sSatisfies the following formula:
Figure FDA0003612457410000018
wherein a, b belongs to {0,1, …, K-1}, m belongs to { - (2Q-2), - (2Q-3), …,2Q-2},
Figure FDA0003612457410000019
fsis the sampling rate, 0 < fs＜fmax，fmaxFor signal frequency upper limit, when the sampling value is interfered by noise, in order to improve the robustness of the method, the sampling rate f 'generated by the sampling rate generator module is fed back'sThe following equation needs to be satisfied:
Figure FDA0003612457410000021
wherein rem { (. 1) } represents a remainder of division of (-) by 1, and ε represents a threshold determined by noise intensity;
step seven, the switch is set by the control signal to be in a sub-sampling channel, the sub-sampling channel samples the obtained shunt signal, a sampling clock CLK2 samples the shunt signal uniformly, and the sampling rate is f'sThe sample values are represented as follows:
Figure FDA0003612457410000022
wherein, c'kIn order to be a complex amplitude of the signal,
Figure FDA0003612457410000023
T′s＝1/f′sis the sampling period, σ [ n']For model matching the sampled values of the error signal, N 'is 0,1, …, N' -1, and then an optimal estimation algorithm of the number of amplitude and frequency components is run, the sampled values x 'N']By the algorithm, an estimated value is obtained
Figure FDA0003612457410000024
Step eight, estimating the parameter
Figure FDA0003612457410000025
And
Figure FDA0003612457410000026
input to a joint estimation parameter algorithm for complex frequency and complex amplitude
Figure FDA0003612457410000027
And
Figure FDA0003612457410000028
and (4) obtaining a joint estimation.
2. The method for measuring the parameters of the non-ideal multi-damping harmonic signal based on the two-channel undersampling as claimed in claim 1, wherein the process of the fourth step is as follows:
step 4.1, the established optimization model is expressed as follows:
Figure FDA0003612457410000029
due to the fact that
Figure FDA00036124574100000210
(3) To convert to:
Figure FDA00036124574100000211
step 4.2, because the input non-ideal med signal x (t) is unknown, that is, the formula (4) cannot be solved, and the sampling value obtained by the formula (2) is used to replace x (t), the formula (4) is converted into the following formula:
Figure FDA00036124574100000212
wherein, the first and the second end of the pipe are connected with each other,
Figure FDA00036124574100000213
number of frequency components K and parameter
Figure FDA00036124574100000214
Are all unknown, then equation (5) is described by:
Figure FDA00036124574100000215
wherein, the first and the second end of the pipe are connected with each other,
Figure FDA0003612457410000031
equation (7) is referred to as the optimization objective function.
3. The method for measuring parameters of a non-ideal multi-damping harmonic signal based on two-channel undersampling as claimed in claim 1 or 2, characterized in that the procedure of the fifth step is as follows:
step 5.1, initialization, first, obtaining the sampling value x [ n ] from the parameter measurement system]N is 0,1, …, N-1, the number of frequency components K' is set to 0, the optimum number of frequency components K is set to 0, and the estimation parameter is set to PbestSetting the optimum value of equation (7) to fbest；
Step 5.2, updating the number K '═ K' +1 of the frequency components;
step 5.3, using the sampled value x [ n ]]And solving the formula (2) by a spectrum estimation method to estimate the parameters
Figure FDA0003612457410000032
As an initial bit for the next step;
step 5.4, the optimization problem of the formula (6) is solved by using the provided optimization algorithm, and P is setoptIs a global optimum position;
step 5.5, update the optimal solution
Figure FDA0003612457410000033
Step 5.6, the result is output, if
Figure FDA0003612457410000034
And (4) outputting to obtain the optimal solution, otherwise, returning to the step 5.3 and repeating the steps.
4. The method for measuring parameters of a non-ideal multi-damping harmonic signal based on two-channel undersampling as claimed in claim 1 or 2, characterized in that the procedure of the step eight is as follows:
step 8.1, estimated normalized frequency due to the periodicity of the trigonometric function
Figure FDA0003612457410000035
And
Figure FDA0003612457410000036
by 2 pi mkNamely:
Figure FDA0003612457410000037
further, the following formula is obtained:
Figure FDA0003612457410000038
wherein, the angle is the main value of the complex angle, and is more than or equal to 0 and less than 2 pi, mkE.g., Z, K e {1,2, …, K }, since
Figure FDA0003612457410000039
fmaxFor the upper frequency limit of the MEDS signal, the following equation is obtained:
Figure FDA00036124574100000310
note book
Figure FDA00036124574100000311
Then (13) is expressed as the following equation:
0≤mk≤(Q-1) (14)
step 8.2, according to (12), from
Figure FDA00036124574100000312
The argument principal value of (a) is obtained as a frequency parameter solution set F and expressed as follows:
Figure FDA0003612457410000041
wherein, fsIs the sampling rate of the main sampling channel, the angle is the amplitude main value of the complex number, the angle is more than or equal to 0 and less than 2 pi, mp∈Z，-(Q-1)≤nk≤Q-1，nkE.g. Z, so the true frequency
Figure FDA0003612457410000042
Meanwhile, the estimated value of the damping factor is obtained as follows:
Figure FDA0003612457410000043
step 8.3, in the sub-sampling channel, according to (12), from
Figure FDA0003612457410000044
The argument principal value of (a) is obtained as a frequency parameter solution set F' and expressed in the following form:
Figure FDA0003612457410000045
wherein, f'sIs the sampling rate of a secondary sampling channel, and the angle () is the amplitude primary value of a complex number (-), and is more than or equal to 0 and less than 2 pi, m'p′∈Z，-(Q′-1)≤n′k≤Q′-1，n′kE.g. Z, thus, the true frequency
Figure FDA0003612457410000046
Meanwhile, the estimated value of the damping factor is obtained as follows:
Figure FDA0003612457410000047
thus, the frequency parameter is obtained as shown in the following formula:
Figure FDA0003612457410000048
and the damping factor parameter is as follows:
Figure FDA0003612457410000049
thereby obtaining a complex frequency skThe estimate is shown below:
Figure FDA00036124574100000410
step 8.4, once the frequency is found
Figure FDA00036124574100000411
Calculating the parameters
Figure FDA00036124574100000412
And
Figure FDA00036124574100000413
first, by selecting
Figure FDA00036124574100000414
And
Figure FDA00036124574100000415
the same elements in (1), assuming other elements are
Figure FDA00036124574100000416
And
Figure FDA00036124574100000417
then, the parameters
Figure FDA00036124574100000418
And
Figure FDA00036124574100000419
respectively from the set
Figure FDA00036124574100000420
And
Figure FDA00036124574100000421
is obtained, therefore, parameter ckEstimated as:
Figure FDA00036124574100000422
CN202110103500.8A 2021-01-26 2021-01-26 Non-ideal multi-damping harmonic signal parameter measurement method based on two-channel undersampling Active CN112964931B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110103500.8A CN112964931B (en) 2021-01-26 2021-01-26 Non-ideal multi-damping harmonic signal parameter measurement method based on two-channel undersampling

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110103500.8A CN112964931B (en) 2021-01-26 2021-01-26 Non-ideal multi-damping harmonic signal parameter measurement method based on two-channel undersampling

Publications (2)

Publication Number Publication Date
CN112964931A CN112964931A (en) 2021-06-15
CN112964931B true CN112964931B (en) 2022-07-15

Family

ID=76272648

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110103500.8A Active CN112964931B (en) 2021-01-26 2021-01-26 Non-ideal multi-damping harmonic signal parameter measurement method based on two-channel undersampling

Country Status (1)

Country Link
CN (1) CN112964931B (en)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115276799B (en) * 2022-07-27 2023-07-11 西安理工大学 Decision threshold self-adaption method for undersampling modulation demodulation in optical imaging communication

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101713795A (en) * 2009-09-09 2010-05-26 中国科学院国家授时中心 Method of digitalized measuring frequency in restriction of sampling rate
KR20150047109A (en) * 2013-10-23 2015-05-04 삼성전자주식회사 Magnetic resonance imaging apparatus and imaging method for magnetic resonance image thereof
EP3301461A1 (en) * 2016-09-28 2018-04-04 Siemens Aktiengesellschaft Method for detection of harmonics of a univariate signal
CN111224672A (en) * 2020-01-16 2020-06-02 哈尔滨工业大学 Multi-harmonic signal undersampling method based on multi-channel time delay

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101713795A (en) * 2009-09-09 2010-05-26 中国科学院国家授时中心 Method of digitalized measuring frequency in restriction of sampling rate
KR20150047109A (en) * 2013-10-23 2015-05-04 삼성전자주식회사 Magnetic resonance imaging apparatus and imaging method for magnetic resonance image thereof
EP3301461A1 (en) * 2016-09-28 2018-04-04 Siemens Aktiengesellschaft Method for detection of harmonics of a univariate signal
CN111224672A (en) * 2020-01-16 2020-06-02 哈尔滨工业大学 Multi-harmonic signal undersampling method based on multi-channel time delay

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
三通道欠采样频率估计;黄山等;《华中科技大学学报(自然科学版)》;20170930(第09期);第11-15页 *

Also Published As

Publication number Publication date
CN112964931A (en) 2021-06-15

Similar Documents

Publication Publication Date Title
CN107102255B (en) Single ADC acquisition channel dynamic characteristic test method
CN107192878A (en) A kind of trend of harmonic detection method of power and device based on compressed sensing
CN109407501B (en) Time interval measuring method based on relevant signal processing
CN106546817B (en) A kind of Frequency Estimation and energy state postulate with feedback function
Zou et al. Blind timing skew estimation using source spectrum sparsity in time-interleaved ADCs
CN111224672B (en) Multi-channel delay-based multi-harmonic signal undersampling method
CN112964931B (en) Non-ideal multi-damping harmonic signal parameter measurement method based on two-channel undersampling
CN105720983A (en) Error estimation method and device for time interleaving analog-digital conversion system
CN109889231A (en) Burst signal lack sampling method based on random demodulation and the limited new fixed rate of interest
Fu et al. Sub-Nyquist sampling of multiple sinusoids
CN112859019A (en) Intra-pulse modulation type parameter extraction system and using method
CN112162152B (en) Sine wave coherent pulse train signal frequency estimation method based on phase straight line fitting
CN105182069A (en) High resolution group quantization phase processing method under pilot frequency architecture
Fu et al. Parameter Measurement of $ M $-Ary PSK Signals With Finite Rate of Innovation
Chen et al. Robust precise time difference estimation based on digital zero-crossing detection algorithm
CN103543331B (en) A kind of method calculating electric signal harmonic wave and m-Acetyl chlorophosphonazo
Djurovic Estimation of the sinusoidal signal frequency based on the marginal median DFT
CN110808929A (en) Real-complex conversion type signal-to-noise ratio estimation algorithm of subtraction strategy
CN114301552B (en) Digital modulation signal testing method and system
CN110944336A (en) Time-frequency spectrum sensing method based on limited new information rate
Ta et al. Fully digital background calibration technique for channel mismatches in TIADCs
Kusuma et al. On the accuracy and resolution of powersum-based sampling methods
CN112595889B (en) under-Nyquist sampling and parameter measuring method for non-ideal multi-exponential decay sinusoidal signal
Huang et al. Sub-Nyquist sampling of multiple exponentially damped sinusoids with feedback structure
CN112883787B (en) Short sample low-frequency sinusoidal signal parameter estimation method based on spectrum matching

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