%MVDR最小方差无失真响应
clear;clc;
%spatially resample
c = 1500;
f = [84:120];
d = c/f(end)/2;
K =7;
M = 2*K+1; %15 阵元个数
P = 2; %信号源个数
u = [86 96];
fs = 2*f(end); %240
snrdB=12;
snr=10^(snrdB/20);
snapshots = M*(fs+1); % 3615
t = [0:snapshots-1]/fs; %0:15.0583
% for trial = 1:5
%!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!%generate source
phi0 = randn(P,length(f))*2; %信号各频率初相位
df = [0 0]; %各源频率差
for kf=1:length(f)
for k=1:P
s(k,:,kf) = snr*exp(j*( 2*pi*( f(kf)+df(k) )*t +(k-1)*0.0+phi0(1,kf)));%源
a(:,k,kf) = exp(j*( 2*pi*( f(kf) )*cos(u(k)*pi/180)*d/c*[-K:K]' ));%流型
end
n(:,:,kf) = exp(j*2*pi*randn(M,snapshots));%噪声
end
for kf=1:length(f)
x(:,:,kf) = a(:,:,kf)*s(:,:,kf)+n(:,:,kf);
end
% clear s;clear a;clear n;
% %!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!%mvdr
% theta = [0:180];mvdr=zeros(1,length(theta));
% for kf = 1:length(f)
% Rf = x(:,:,kf)*x(:,:,kf)'/snapshots;
% for k = 1:length(theta)
% Dtheta=exp(-j*2*pi*f(kf)*cos(theta(k)*pi/180)*d/c*[0:M-1].');
% R=Dtheta'*inv(Rf)*Dtheta;
% mvdr(k) = mvdr(k)+1/R;
% end
% end
% plot(theta,10*log10(abs(mvdr)));grid on;
% %!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!%mvdr
% sx = zeros(M,snapshots);
% for kf=1:length(f)
% sx = a(:,:,kf)*s(:,:,kf)+sx;
% end
% clear s;clear a;clear n;
% x = sx+1*exp(j*2*pi*randn(size(sx)));clear sx;
%!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!%合成
temp = zeros(M,snapshots);
for kf=1:length(f)
temp = x(:,:,kf)+temp;
end
x = temp; %size(x)=15*3615
clear temp;
clear s;clear a;clear n;
%!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!%
nfft = fs+1; %241
fL = f(1);fH = f(end);
dtheta = 0.1;
theta=0:dtheta:180;
f_L=(nfft-1)/fs*fL+1; %85
f_H=(nfft-1)/fs*fH+1; %121
Bn=f_H-f_L+1; %37
[M Lx]=size(x); %15*3615
N = fix(Lx/nfft); %N=15
%!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!%互谱密度矩阵
R_csdm_f=zeros(M,M*Bn);%MxM matrix for each frequency bin 15*555
for ii=1:N
seg=x(:,(ii-1)*(nfft)+1:(ii)*(nfft)); %15*241
fftseg=fft(seg.',nfft); %241*15
for i=f_L:f_H
KK(:,(i-f_L)*M+1:(i-f_L+1)*M)=fftseg(i,:).'*conj(fftseg(i,:)); %size(KK)=15*555
end
R_csdm_f=R_csdm_f+KK;
end
R_csdm_f=1/N*R_csdm_f;
clear i ii;
%!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!%mvdr
op_mvdr=zeros(1,length(theta));
for i=1:length(theta)
for ii=1:Bn
Dtheta=exp(j*2*pi*(ii-1+f_L-1)*fs/(nfft-1)*cos(theta(i)*pi/180)*d/c*[0:M-1].');
R2=Dtheta'*inv(R_csdm_f(:,(ii-1)*M+1:ii*M))*Dtheta;
op_mvdr(i)=1/R2+op_mvdr(i);
end
end
hold on;
plot(theta,10*log10(abs(op_mvdr)))
% end
%!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!%srmv
f0 = f(1);%聚焦频率 84
R_srmv = zeros(M,M);
for k = f_L:f_H
fi = k-1; %84:
fai=min(pi,(fi/f0)*pi);
Ky=round(K*f0/f0) %7
T = zeros(2*Ky+1,M); %15*15
for n=-Ky:Ky
for m=-K:K
if ((f0/fi)*n-m==0)
T(Ky+1+n,K+1+m) = fai/pi;
else
T(Ky+1+n,K+1+m)=sin(fai*((f0/fi)*n-m))/(((f0/fi)*n-m)*pi);
end
end
end
R_srmv = R_srmv+T*R_csdm_f( :,(k-f_L)*M+1:(k-f_L+1)*M )*T';
end
for k=1:length(theta)
Dtheta=exp(j*2*pi*f0*cos(theta(k)*pi/180)*d/c*[-K:K].');
srmv(k)=1/(Dtheta'*inv(R_srmv)*Dtheta);
end
hold on;plot(theta,10*log10(abs(srmv)),'r')
plotxline(u,'m','--')
%!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!%srmv
54225024MVDR.rar_matlab例程_matlab_
版权申诉
5星 · 超过95%的资源 76 浏览量
2021-08-10
02:05:33
上传
评论 2
收藏 2KB RAR 举报
pudn01
- 粉丝: 40
- 资源: 4万+
最新资源
- 2023-04-06-项目笔记 - 第一百五十四阶段 - 4.4.2.152全局变量的作用域-152 -2024.06.04
- 松哥解协议松哥解协议松哥解协议松哥解协议松哥解协议
- 618节日618节日618节日
- tensorflow-gpu-2.9.1-cp37-cp37m-win-amd64.whl
- tensorflow-gpu-2.9.0-cp37-cp37m-win-amd64.whl
- tensorflow-gpu-2.9.0-cp39-cp39-win-amd64.whl
- lcd daimalcd daima
- 电影领域-推荐算法-个性化内容-观影决策-电影推荐小程序.zip
- 电气控制PLC考试题库
- 如何使用MATLAB简介
资源上传下载、课程学习等过程中有任何疑问或建议,欢迎提出宝贵意见哦~我们会及时处理!
点击此处反馈
- 1
- 2
前往页