白噪声下真实正弦波的精确频率估计研究(Matlab代码实现)

2023-12-05

???????????????? 欢迎来到本博客 ❤️❤️????????

????博主优势: ???????????? 博客内容尽量做到思维缜密,逻辑清晰,为了方便读者。

⛳️ 座右铭: 行百里者,半于九十。

???????????? 本文目录如下: ????????????

目录

????1 概述

????2 运行结果

????3 参考文献

????4 Matlab代码及数据


????1 概述

在白噪声下,真实正弦波的精确频率估计是一个重要的研究课题,因为白噪声会对信号的频率进行干扰,使得频率估计变得困难。下面是一些常用的方法来进行真实正弦波的精确频率估计研究:

1. 傅里叶变换:通过对信号进行傅里叶变换,可以得到信号的频谱信息,从而可以估计出频率。然而,在白噪声下,信号的频谱可能会被噪声掩盖,导致频率估计的不准确性。

2. 自相关函数:自相关函数可以用来评估信号中重复出现的模式,从而可以用来估计信号的频率。然而,在白噪声下,自相关函数的峰值可能会被噪声掩盖,导致频率估计的不准确性。

3. 峰值检测:通过寻找信号的峰值来进行频率估计。然而,在白噪声下,峰值可能会被噪声掩盖,导致频率估计的不准确性。

4. 希尔伯特变换:希尔伯特变换可以将信号从时域转换到频域,从而可以得到信号的频率信息。然而,在白噪声下,希尔伯特变换可能会受到噪声的影响,导致频率估计的不准确性。

综上所述,白噪声下真实正弦波的精确频率估计是一个复杂的问题,需要综合运用多种方法来进行研究和解决。在实际应用中,可以根据具体的情况选择合适的方法来进行频率估计。

???? 2 运行结果

可视化代码:

figure()
semilogy(SNRdB,CRLB/Fs*N,'--ks',SNRdB,RMSEBin(:,1),'-mo',SNRdB,RMSEBin(:,2),'-r*',...
SNRdB,RMSEBin(:,3),'-gd',SNRdB,RMSEBin(:,4),'-c.',SNRdB,RMSEBin(:,5),'-b^',...
SNRdB,RMSEBin(:,6),'-bv','LineWidth',1.2)
legend('CRLB','Proposed, second','Proposed, third','PSVD [5]','Djukanovic [9]','Ye [12],\it{T}=2','Ye [12],\it{T}=4')
xlabel('SNR(dB)')
ylabel('RMSE(in DFT bins)')
axis([-10,60,5*10^-5,20])
grid on
saveas(gcf,'..\results\Fig1.png')

%% Fig 2
Fs = 1; % Sampling frequency
N = [16 64 256 1024 4096]; % Number of samples
A = 1; % Half of Amplitude
SimTimes = 2000; % Number of simulation
SNRdB = 20; % SNR
RMSE = zeros(length(N),6); % RMSE
CRLB = sqrt(3*Fs^2./(pi^2*N.*(N.^2-1)*10.^(SNRdB/10))).'; % CRLB

for ii = 1:length(N)
ErrF = zeros(SimTimes,6); % Error of frequency estimation
for jj = 1:SimTimes
F = Fs/N(ii)*1.2; % Frequency
P = rand*2*pi; % Initial Phase
n = (0:N(ii)-1).';
s = 2*A*cos(2*pi*F/Fs*n+P); % real sinusoid
Ps = mean(s.^2); % Signal power
Pn = Ps/(10^(SNRdB/10)); % Noise Power
w = sqrt(Pn)*randn(N(ii),1);
s = w+s; % AWGN channel

EstF = ProposedInitial(s,Fs); % Initial estimate of proposed algorithm
ErrF(jj,1) = abs(EstF-F);
if ErrF(jj,1)>Fs/2
ErrF(jj,1) = Fs-ErrF(jj,1);
end

EstF = ProposedFine(s,Fs); % Fine estimate of proposed algorithm
ErrF(jj,2) = abs(EstF-F);
if ErrF(jj,2)>Fs/2
ErrF(jj,2) = Fs-ErrF(jj,2);
end

EstF = PSVD(s,Fs);
ErrF(jj,3) = abs(EstF-F);
if ErrF(jj,3)>Fs/2
ErrF(jj,3) = Fs-ErrF(jj,3);
end

EstF = Djukanovic(s,Fs);
ErrF(jj,4) = abs(EstF-F);
if ErrF(jj,4)>Fs/2
ErrF(jj,4) = Fs-ErrF(jj,4);
end

EstF = Ye(s,Fs,2); % Ye algorithm with T=2
ErrF(jj,5) = abs(EstF-F);
if ErrF(jj,5)>Fs/2
ErrF(jj,5) = Fs-ErrF(jj,5);
end

EstF = Ye(s,Fs,4); % Ye algorithm with T=4
ErrF(jj,6) = abs(EstF-F);
if ErrF(jj,6)>Fs/2
ErrF(jj,6) = Fs-ErrF(jj,6);
end


end
RMSE(ii,:) = sqrt(mean(ErrF.^2)); % RMSE in DFT bins
end
RMSEBin = RMSE/Fs.*repmat(N.',1,6); % RMSE in DFT bins
figure()
semilogy(log2(N),CRLB/Fs.*N.','--ks',log2(N),RMSEBin(:,1),'-mo',log2(N),RMSEBin(:,2),'-r*',...
log2(N),RMSEBin(:,3),'-gd',log2(N),RMSEBin(:,4),'-c.',log2(N),RMSEBin(:,5),'-b^',...
log2(N),RMSEBin(:,6),'-bv','LineWidth',1.2)
legend('CRLB','Proposed, second','Proposed, third','PSVD [5]','Djukanovic [9]','Ye [12],\it{T}=2','Ye [12],\it{T}=4')
xlabel('log_2\it{N}')
ylabel('RMSE(in DFT bins)')
grid on
saveas(gcf,'..\results\Fig2.png')

%% Fig 3
Fs = 1; % Sampling frequency
N = 1024; % Number of samples
A = 1; % Half of Amplitude
SimTimes = 1000; % Number of simulation
SNRdB = 20; % SNR
x = [0.55:0.25:2.05,3,6:N/2/9:N/2-6,N/2-3,N/2-2.05:0.25:N/2-0.55]; % Frequency in DFT bins
F = x*Fs/N; % Frequency in Hz
RMSE = zeros(length(F),6); % RMSE
CRLB = sqrt(3*Fs^2./(pi^2*N.*(N.^2-1)*10.^(SNRdB/10)))*ones(length(F),1); % CRLB

for ii = 1:length(F)
ErrF = zeros(SimTimes,6);
CRLBInst = zeros(SimTimes,1);
for jj = 1:SimTimes
P = rand*2*pi; % Initial Phase
n = (0:N-1).';
s = 2*A*cos(2*pi*F(ii)/Fs*n+P); % real sinusoid
Ps = mean(s.^2); % Signal power
Pn = Ps/(10^(SNRdB/10)); % Noise Power
w = sqrt(Pn)*randn(N,1);
s = w+s; % AWGN channel

EstF = ProposedInitial(s,Fs); % Initial estimate of proposed algorithm
ErrF(jj,1) = abs(EstF-F(ii));
if ErrF(jj,1)>Fs/2
ErrF(jj,1) = Fs-ErrF(jj,1);
end

EstF = ProposedFine(s,Fs); % Fine estimate of proposed algorithm
ErrF(jj,2) = abs(EstF-F(ii));
if ErrF(jj,2)>Fs/2
ErrF(jj,2) = Fs-ErrF(jj,2);
end

EstF = PSVD(s,Fs);
ErrF(jj,3) = abs(EstF-F(ii));
if ErrF(jj,3)>Fs/2
ErrF(jj,3) = Fs-ErrF(jj,3);
end

EstF = Djukanovic(s,Fs);
ErrF(jj,4) = abs(EstF-F(ii));
if ErrF(jj,4)>Fs/2
ErrF(jj,4) = Fs-ErrF(jj,4);
end

EstF = Ye(s,Fs,2); % Ye algorithm with T=2
ErrF(jj,5) = abs(EstF-F(ii));
if ErrF(jj,5)>Fs/2
ErrF(jj,5) = Fs-ErrF(jj,5);
end

EstF = Ye(s,Fs,4); % Ye algorithm with T=4
ErrF(jj,6) = abs(EstF-F(ii));
if ErrF(jj,6)>Fs/2
ErrF(jj,6) = Fs-ErrF(jj,6);
end

end
RMSE(ii,:) = sqrt(mean(ErrF.^2,1));
end
RMSEBin = RMSE/Fs*N; % RMSE in DFT bins
figure()
semilogy(x,CRLB/Fs*N,'--ks',x,RMSEBin(:,1),'-mo',x,RMSEBin(:,2),'-r*',...
x,RMSEBin(:,3),'-gd',x,RMSEBin(:,4),'-c.',x,RMSEBin(:,5),'-b^',...
x,RMSEBin(:,6),'-bv','LineWidth',1.2)
legend('CRLB','Proposed, second','Proposed, third','PSVD [5]','Djukanovic [9]','Ye [12],\it{T}=2','Ye [12],\it{T}=4')
xlabel('\it{Nf/fs}')
ylabel('RMSE(Hz)')
axis([0,512,1*10^-3,1])
grid on
saveas(gcf,'..\results\Fig3.png')

????3 参考文献

文章中一些内容引自网络,会注明出处或引用为参考文献,难免有未尽之处,如有不妥,请随时联系删除。

[1] Ying X , Shu-Xiong L , Jun-Xun Y .Frequency estimation of single complex sinusoid in white Gaussian noise[J].Journal of China Institute of Communications, 2002.DOI:10.1002/mop.10502.

[2] Liu S , Xu W , Wang Z .An Accurate Frequency Estimator of Single Sinusoid in Gaussian Colored Noise.[C]//International Conference on Innovative Computing.IEEE, 2006.DOI:10.1109/ICICIC.2006.227.

[4]丁志中.白噪声背景下正弦信号频率和功率估计的有效算法[J].合肥工业大学学报:自然科学版, 1995, 18(4):6.DOI:CNKI:SUN:HEFE.0.1995-04-002.

[5]吴珊珊胡国兵丁宁汤滟.正弦波频率估计结果的可靠性评估算法研究[J].现代雷达, 2015, 37(3):31-35.

???? 4 Matlab代码及数据

本文内容由网友自发贡献,版权归原作者所有,本站不承担相应法律责任。如您发现有涉嫌抄袭侵权的内容,请联系:hwhale#tublm.com(使用前将#替换为@)

白噪声下真实正弦波的精确频率估计研究(Matlab代码实现) 的相关文章

  • 将 Matlab 数组移植到 C/C++

    我正在将 matlab 程序移植到 C C 我有几个问题 但最重要的问题之一是 Matlab 将任何维度的数组都视为相同 假设我们有一个这样的函数 function result f A B C result A 2 B C A B and
  • 在 MATLAB 中绘图后恢复轴

    从文本文件绘制多种方法的输出后 未显示轴的右侧和上侧 我需要拥有它们并将它们加粗 就像当前的轴一样 绘制的数据来自存储每种方法数据的文件 每个数据文件都是一个 256x2 文件 包含 0 1 之间的值 第一列是精度 第二列是召回率 figu
  • 通过 cuFFT 进行逆 FFT 缩放

    每当我使用 cuFFT 绘制程序获得的值并将结果与 Matlab 的结果进行比较时 我都会得到相同形状的图形 并且最大值和最小值位于相同的点 然而 cuFFT 得到的值比 Matlab 得到的值大得多 Matlab代码是 fs 1000 s
  • 如何在 MATLAB 中将矩阵元素除以列总和?

    有没有一种简单的方法可以将每个矩阵元素除以列和 例如 input 1 4 4 10 output 1 5 4 14 4 5 10 14 以下是执行此操作的不同方法的列表 使用bsxfun https www mathworks com he
  • matlab 中的动画绘图

    我正在尝试创建一个三角形的动画图 最终结果应该是十个三角形 后面跟着两个更大的三角形 后面跟着一条直线 使用matlab文档 https de mathworks com help matlab ref drawnow html 我最终得到
  • 保存符号方程以供以后使用?

    From here http www mathworks com help releases R2011a toolbox symbolic brvfu8o 1 html brvfxem 1 我正在尝试求解这样的符号方程组 syms x y
  • 考虑预分配速度[重复]

    这个问题在这里已经有答案了 我正在做以下事情 for i 1 m index 0 for j 1 n index index values i j 2 j 1 if j 1 symbol chip chip values index 1 e
  • 在矩阵中找到叉的最快方法

    定义 A i j 1 是十字的中点 如果元素A i 1 j 1A i 1 j 1A i j 1 1A i j 1 1 这些元素和中点一起形成矩阵 A 中的十字 其中 A 至少是一个 3 3 矩阵 并且i j 0 假设上图是 8 8 矩阵 A
  • 氡变换线检测

    我正在尝试检测灰度图像中的线条 为此 我在 MATLAB 中使用 Radon 变换 我的 m 文件的示例如下所示 我可以使用此代码检测多行 我还使用线条的移位和旋转属性来绘制线条 但是 我不明白在获取rho和theta值后如何获取检测线的起
  • MATLAB:图像角坐标和引用元胞数组

    我在比较不同元胞数组中的元素时遇到一些问题 这个问题的背景是我正在使用bwboundariesMATLAB 中的函数可追踪图像的轮廓 该图像是结构横截面 我试图找出整个部分是否具有连续性 即 只有一个轮廓由bwboundaries命令 完成
  • 在 Python 上显示 Matlab mat 文件中的图像

    我目前正在尝试显示从此下载的 Mat 文件中的图像site http www rctn org bruno sparsenet 这是一个 mat 文件 所以我尝试使用 scipy io loadmat 函数加载它 但我似乎无法绘制图像 我究
  • matlab中优先级队列的实现方法

    matlab中有没有提供minpriorityqueue功能的库 import java util PriorityQueue import java util public class MyQueue Comparator
  • 在Matlab图例中使用Latex?

    我的 matlab 不接受我的 Latex 例如 如果我使用legend b 6 rightarrow b 7 它没有向我显示箭头 我该如何解决这个问题 尝试使用 Latex 解释器 例如 legend b 6 rightarrow b 7
  • 使用网络计算机进行 Matlab 并行处理

    我熟悉matlabpool and parfor用法 但我仍然需要加快计算速度 我的 1GB 网络中有一台功能更强大的计算机 两台计算机都有 R2010b 并且具有相同的代码和路径 使用两台计算机进行并行计算的最简单方法是什么 我今天使用的
  • 使用 scipy.io 将 python pandas dataframe 转换为 matlab 结构

    我正在尝试使用 scipy io 将 pandas 数据帧保存到 matlab mat 文件 我有以下内容 array1 np array 1 2 3 array2 np array a b c array3 np array 1 01 2
  • 使用正常数据直方图与直接公式进行熵估计(matlab)

    假设我们已经绘制了n 10000标准正态分布的样本 现在我想使用直方图计算其熵来计算概率 1 计算概率 例如使用matlab p x hist samples binnumbers area x 2 x 1 sum p p p area b
  • Matlab的uicontrol在Octave中的实现?

    我正在尝试在 Octave 中运行我们实验室中使用的图形程序的 m Matlab 代码 Octave 告诉我代码中使用的函数 uicontrol 没有定义 经过一番搜索 我发现 JHandles 包有一个 uicontrol GUI 功能的
  • 可以避免迭代元胞数组时的“s{1} 烦恼”吗?

    The s 1 标题的 烦恼 指的是下面的 for 块中的第一行 for s some cell array s s 1 unpeel the enclosing cell do stuff with s end This s s 1 业务
  • 静态时序数据的数据库解决方案

    我们拥有一个庞大且不断增长的实验数据集 该数据集取自约 30 000 名受试者 对于每个主题 都有多个数据记录 在每个记录中 收集了多个生理数据时间序列 每个时间序列约 90 秒长 并以 250Hz 采样 我应该注意到 时间序列的任何给定实
  • MATLAB:MEX 矩阵除法给出的结果与 m 文件不同

    我使用 MATLAB 的编码器工具创建了矩阵指数函数的 MEX 版本 以在另一组函数中使用 问题是 MEX 版本给出的结果与原始 m 文件不同 经过调试 我认为这是因为MEX文件和m文件没有做相同的矩阵除法 或者 MEX 文件首先就有问题

随机推荐