clear,clc
m=8; % sensors
n=2; % sources
theta=[-20 0]; % in angle
d=1/2; % 1/2 lambada
N=100; % samples
L=100; % resolution in [-90' 90']
meann=0; % mean of noise
varn=1; % variance of noise
SNR=10;
INR=20;
rvar1=sqrt(varn) * 10^(SNR/20); % variance of signal
rvar2=sqrt(varn) * 10^(INR/20); % variance of interference
% generate the source signals
s=[rvar1*exp(j*2*pi*(50*0.001*[0:N-1]+rand))
rvar2*exp(j*2*pi*(100*0.001*[0:N-1]+rand))];
% generate the A matrix
A=exp(-j*2*pi*d*[0:m-1].'*sin(theta*pi/180));
% generate the noise component
e=sqrt(varn/2)*(randn(m,N)+j*randn(m,N));
% generate the ULA data
Y=A*s+e;
% compute R
R=Y*Y'/N;
% inverse of R
IR=inv(R);
% beamform using the MVDR method
a1=A(:,1);
w=IR*a1/real(a1'*IR*a1);
beam=zeros(1,L);
for i = 1 : L
a=exp(-j*2*pi*d*[0:m-1].'*sin(-pi/2 + pi*(i-1)/L));
beam(i)=20*log10(abs(w'*a));
end
angle=-90:180/L:(90-180/L);
plot(angle,beam,'b');
hold on;
%*****************************************************************
m=8;
n=2;
theta=[-20 0];
d=1/2;
N=100;
L=100;
meann=0;
varn=1;
SNR=10;
INR=30;
rvar1=sqrt(varn) * 10^(SNR/20);
rvar2=sqrt(varn) * 10^(INR/20);
s=[rvar1*exp(j*2*pi*(50*0.001*[0:N-1]+rand))
rvar2*exp(j*2*pi*(100*0.001*[0:N-1]+rand))];
A=exp(-j*2*pi*d*[0:m-1].'*sin(theta*pi/180));
e=sqrt(varn/2)*(randn(m,N)+j*randn(m,N));
Y=A*s+e;
R=Y*Y'/N;
IR=inv(R);
% beamform using the MVDR method
a1=A(:,1);
w=IR*a1/real(a1'*IR*a1);
beam=zeros(1,L);
for i = 1 : L
a=exp(-j*2*pi*d*[0:m-1].'*sin(-pi/2 + pi*(i-1)/L));
beam(i)=20*log10(abs(w'*a));
end
angle=-90:180/L:(90-180/L);
plot(angle,beam,'r');
%****************************************************************
m=8;
n=2;
theta=[-20 0];
d=1/2;
N=100;
L=100;
meann=0;
varn=1;
SNR=10;
INR=40;
rvar1=sqrt(varn) * 10^(SNR/20);
rvar2=sqrt(varn) * 10^(INR/20);
s=[rvar1*exp(j*2*pi*(50*0.001*[0:N-1]+rand))
rvar2*exp(j*2*pi*(100*0.001*[0:N-1]+rand))];
A=exp(-j*2*pi*d*[0:m-1].'*sin(theta*pi/180));
e=sqrt(varn/2)*(randn(m,N)+j*randn(m,N));
Y=A*s+e;
R=Y*Y'/N;
IR=inv(R);
% beamform using the MVDR method
a1=A(:,1);
w=IR*a1/real(a1'*IR*a1);
beam=zeros(1,L);
for i = 1 : L
a=exp(-j*2*pi*d*[0:m-1].'*sin(-pi/2 + pi*(i-1)/L));
beam(i)=20*log10(abs(w'*a));
end
angle=-90:180/L:(90-180/L);
plot(angle,beam,'g');
xlabel('angle');
ylabel('幅度响应/dB');
gtext({'\bf\leftarrowINR=20';'\bf\leftarrowINR=30';'\bf\leftarrowINR=40'})