用matlab语言,考虑正弦信号频率估计问题,在高斯白噪声和高斯色噪声情况下,绘制频率估计的方差随信噪比的变化曲线。其中高斯色噪声由AR1模型产生,模型参数为0.8,运用蒙特卡洛方法,频率真值为0.25,运用最大似然法进行频率估计
时间: 2024-06-11 13:10:54 浏览: 4
以下是代码实现:
%% 参数设置
fs = 1000; % 采样频率
T = 1; % 信号时长
N = T*fs; % 采样点数
f_true = 0.25; % 真实频率
SNR_db = -10:2:40; % 信噪比范围(dB)
num_trials = 100; % 蒙特卡洛试验次数
AR_coef = 0.8; % AR1模型参数
%% 信号生成
t = (0:N-1)/fs;
x = sin(2*pi*f_true*t);
%% 蒙特卡洛试验
variance_ml = zeros(size(SNR_db)); % 存储最大似然法频率估计的方差
variance_crlb = zeros(size(SNR_db)); % 存储CRLB计算结果的方差
for i = 1:length(SNR_db)
snr = 10^(SNR_db(i)/10); % 信噪比
noise_power = norm(x)^2/(N*snr); % 噪声功率
for j = 1:num_trials
% 生成噪声
white_noise = sqrt(noise_power)*randn(1,N); % 高斯白噪声
ar1_noise = filter(1,[1 -AR_coef],white_noise); % 高斯色噪声
% 信号加噪声
y = x + ar1_noise;
% 最大似然法估计频率
[f_ml,~] = freq_est_ml(y,fs);
% 计算CRLB
crlb = freq_est_crlb(y,fs,f_true);
% 存储方差
variance_ml(i) = variance_ml(i) + (f_ml - f_true)^2;
variance_crlb(i) = variance_crlb(i) + crlb;
end
variance_ml(i) = variance_ml(i)/num_trials;
variance_crlb(i) = variance_crlb(i)/num_trials;
end
%% 绘图
figure;
semilogy(SNR_db,variance_ml,'o-',SNR_db,variance_crlb,'x-');
grid on;
xlabel('SNR (dB)');
ylabel('方差');
legend('最大似然法','CRLB');
title(['AR1模型参数:',num2str(AR_coef)]);
%% 最大似然法频率估计函数
function [f_ml,likelihood] = freq_est_ml(x,fs)
x = x(:); % 确保输入为列向量
N = length(x);
X = fft(x);
L = abs(X(2:floor(N/2))).^2;
[~,k] = max(L);
f_ml = (k-1)/N*fs; % 最大似然法估计的频率
if nargout > 1 % 如果需要输出似然函数值
likelihood = L(k);
end
end
%% CRLB计算函数
function crlb = freq_est_crlb(x,fs,f_true)
x = x(:); % 确保输入为列向量
N = length(x);
X = fft(x);
L = abs(X(2:floor(N/2))).^2;
k = round(f_true*N/fs)+1;
dL_df = (L(k+1)-L(k-1))/(2*N/fs);
crlb = 1/(2*pi^2*fs^2*dL_df^2); % CRLB计算公式
end
相关推荐
![rar](https://img-home.csdnimg.cn/images/20210720083606.png)
![rar](https://img-home.csdnimg.cn/images/20210720083606.png)
![](https://csdnimg.cn/download_wenku/file_type_ask_c1.png)
![](https://csdnimg.cn/download_wenku/file_type_ask_c1.png)
![](https://csdnimg.cn/download_wenku/file_type_ask_c1.png)
![](https://csdnimg.cn/download_wenku/file_type_ask_c1.png)
![](https://csdnimg.cn/download_wenku/file_type_ask_c1.png)
![](https://csdnimg.cn/download_wenku/file_type_ask_c1.png)
![](https://csdnimg.cn/download_wenku/file_type_ask_c1.png)
![](https://csdnimg.cn/download_wenku/file_type_ask_c1.png)
![](https://csdnimg.cn/download_wenku/file_type_ask_c1.png)
![](https://csdnimg.cn/download_wenku/file_type_ask_c1.png)
![](https://csdnimg.cn/download_wenku/file_type_ask_c1.png)
![](https://csdnimg.cn/download_wenku/file_type_ask_c1.png)
![pdf](https://img-home.csdnimg.cn/images/20210720083512.png)
![rar](https://img-home.csdnimg.cn/images/20210720083606.png)