Capon 波束形成matlab仿真(附源代码)

合集下载

CAPON波束形成_Matlab程序

C A P O N波束形成_M a t l a b程序(总3页)--本页仅作为文档封面，使用时请直接删除即可----内页可以根据需求调整合适字体及大小--CAPON波束形成器仿真1．实验原理波束形成就是从传感器阵列重构源信号。

（1）、通过增加期望信源的贡献来实现；（2）、通过抑制掉干扰源来实现。

经典的波束形成需要观测方向(期望信源的方向)的知识。

盲波束形成试图在没有期望信源方向信息的情况下进行信源的恢复。

波束形成技术的基本思想是：通过将各阵元输出进行加权求和，在一时间内将天线阵列波束“导向”到一个方向上，对期望信号得到最大输出功率的导向位置即给出DOA估计。

虽然阵列天线的方向图是全方向的，但阵列的输出经过加权求和后，却可以被调整到阵列接收的方向增益聚集在一个方向上，相当于形成了一个”波束”。

这就是波束形成的物理意义所在。

在智能天线中，波束形成是关键技术之一，是提高信噪比、增加用户容量的保证，能够成倍地提高通信系统的容量，有效地抑制各种干扰，并改善通信质量。

波束形成器的最佳权向量w取决于阵列方向向量)(ak，而在移动通信里用户的方向向量一般未知，需要估计(称之为DOA估计)。

因此，在计算波束形成的最佳权向量之前，必须在已知阵列几何结构的前提下先估计期望信号的波达方向。

Capon 波束形成器求解的优化问题可表述为w arg min P(w)θ=其约束条件为 1)(a w H =θCapon 波束形成器在使噪声和干扰所贡献的功率为最小的同时，保持了期望信号的功率不变。

因此，它可以看作是一个尖锐的空间带通滤波器。

最优加权向量w 可以利用Lagrange 乘子法求解，其结果为 )(a R ˆ)(a )(a R ˆw 1H 1CAP θθθ--= 当μ不取常数，而取作)(a R ˆ)(a 11H θθμ-=时，最佳权向量就转变成Capon 波束形成器的权向量。

空间谱为)(a R ˆ)(a 1)(P 1-H CAP θθθ=2．变量定义 M ：均匀线阵列数目P ：信号源个数nn ：快拍数angle1、angle2、angle3：信号来波角度u ：复高斯噪声Ps ：信号能量refp ：信噪比（实值）X：接收信号Rxx：接收信号的相关矩阵doa：波达方向估计3．仿真结果采用上述算法进行仿真，结果如图所示。

波束形成Matlab程序

1?均匀线阵方向图%8阵元均匀线阵方向图，来波方向为clc; clear all; close all;0度imag=sqrt(_1);element_num=8;% 阵元数为8d_lamda=1/2;%阵元间距d与波长lamda的关系theta=li nspace(-pi/2,pi/2,200);theta0=0;% 来波方向w=exp(imag*2*pi*dl_lamda*sin(theta0)*[0:eleme nt_nu m-1]');for j=1:le ngth(theta)a=exp(imag*2*pi*dd_l amda*si n(theta(j))*[0:eleme nt_nu m-1]'); p(j)=w'*a;end figure;plot(theta,abs(p)),grid on xlabel('theta/radia n')ylabel('amplitude')title('8 阵元均匀线阵方向图')°2 8阵元均匀线阵方向图7654321-15 -1 -0 5 0 06thetaradian1 158当来波方向为45度时，仿真图如下8阵元均匀线阵方向图如下，来波方向为0 度，20log (dB)8阵元均苛銭阵方向图来波方向为0度随着阵元数的增加，波束宽度变窄，分辨力提高：仿真图如下Q d pE =ro 二2. 波束宽度与波达方向及阵元数的关系clcclear allclose allima=sqrt(-1);element_num1=16; %阵元数element_num2=128;element_num3=1024;lamda=0.03; %波长为0.03 米d=1/2*lamda; % 阵元间距与波长的关系theta=0:0.5:90; forj=1:length(theta);fai(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lamda/(element_n um1*d));psi(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lamda/(element_n um2*d));beta(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lamda/(element_ num3*d));endfigure; plot(theta,fai,'r',theta,psi,'b',theta,beta,'g'),grid on xlabel('theta'); ylabel('Width in radians')title(' 波束宽度与波达方向及阵元数的关系') 仿真图如下：。

波束赋形仿真代码

以下是一个简单的波束赋形仿真代码，使用Python编写：python复制代码import numpy as npimport matplotlib.pyplot as plt# 设置参数num_sensors = 10# 传感器数量num_signals = 5# 信号数量num_time_steps = 100# 时间步数# 生成信号signals = np.random.randn(num_signals, num_time_steps)# 生成传感器位置sensor_positions = np.random.rand(num_sensors, 2) * 10# 在10x10的区域内随机生成位置# 生成波束赋形器参数weights = np.random.randn(num_sensors, num_signals)# 进行波束赋形output = np.dot(sensor_positions, weights.T)# 绘制结果plt.figure()plt.subplot(2, 1, 1)plt.imshow(sensor_positions[:, 0], aspect='auto') # 绘制传感器位置图plt.title('Sensor Positions')plt.xlabel('x')plt.ylabel('y')plt.subplot(2, 1, 2)plt.plot(output) # 绘制波束赋形结果图plt.title('Beamforming Output')plt.xlabel('Time Step')plt.ylabel('Amplitude')plt.show()这个代码使用随机生成信号和传感器位置，并使用随机生成的权重进行波束赋形。

输出结果是一个二维图像，其中第一个子图显示了传感器位置，第二个子图显示了波束赋形结果。

自适应波束形成与Matlab程序代码注解

gridon;
xlabel('theta/radian')
ylabel('amplitude/dB')
title([num2str(element_num)'阵元均匀线阵方向图','来波方向为'num2str(theta0*180/pi)'度']);
axis([-1.5 1.5 -50 0]);
(2)仿真结果
A.来波方向为0°
不归一化
归一化
B.来波方向为45°
不归一化
归一化
C.随着阵元数的增加，波束宽度变窄，分辨力提高，仿真图如下：
非归一化
归一化
不归一化
归一化
2.波束宽度与波达方向及阵元数的关系
（1）matlab 程序
clc;
clearall;
closeall;
imag=sqrt(-1);
element_num1=16;
psi(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lambda/(element_num2*d));
beta(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lambda/(element_num3*d));
end
gridon;
xlabel('theta/radian')
ylabel('amplitude/dB')
title([num2str(element_num)'阵元均匀线阵方向图','来波方向为'num2str(theta0*180/pi)'度']);

Matlab波束形成程序

波束形成与智能天线1.均匀线阵方向图%8阵元均匀线阵方向图，来波方向为0度clc;clear all;close all;imag=sqrt(-1);element_num=8;%阵元数为8d_lamda=1/2;%阵元间距d与波长lamda的关系theta=linspace(-pi/2,pi/2,200);theta0=0;%来波方向w=exp(imag*2*pi*d_lamda*sin(theta0)*[0:element_num-1]');for j=1:length(theta)a=exp(imag*2*pi*d_lamda*sin(theta(j))*[0:element_num-1]');p(j)=w'*a;endfigure;plot(theta,abs(p)),grid onxlabel('theta/radian')ylabel('amplitude')title('8阵元均匀线阵方向图')当来波方向为45度时，仿真图如下：8阵元均匀线阵方向图如下，来波方向为0度，20log（dB）随着阵元数的增加，波束宽度变窄，分辨力提高：仿真图如下：2.波束宽度与波达方向及阵元数的关系clcclear allclose allima=sqrt(-1);element_num1=16; %阵元数element_num2=128;element_num3=1024;lamda=0.03; %波长为0.03米d=1/2*lamda; %阵元间距与波长的关系theta=0:0.5:90;for j=1:length(theta);fai(j)=theta(j)*pi/-asin(sin(theta(j)*pi/)-lamda/(element_num1*d));psi(j)=theta(j)*pi/-asin(sin(theta(j)*pi/)-lamda/(element_num2*d));beta(j)=theta(j)*pi/-asin(sin(theta(j)*pi/)-lamda/(element_num3*d)); endfigure;plot(theta,fai,'r',theta,psi,'b',theta,beta,'g'),grid onxlabel('theta');ylabel('Width in radians')title('波束宽度与波达方向及阵元数的关系')仿真图如下：3. 当阵元间距时，会出现栅瓣，导致空间模糊。

波束形成 Matlab程序

1.均匀线阵方向图%8阵元均匀线阵方向图，来波方向为0度clc;clear all;close all;imag=sqrt(-1);element_num=8;%阵元数为8d_lamda=1/2;%阵元间距d与波长lamda的关系theta=linspace(-pi/2,pi/2,200);theta0=0;%来波方向w=exp(imag*2*pi*d_lamda*sin(theta0)*[0:element_num-1]');for j=1:length(theta)a=exp(imag*2*pi*d_lamda*sin(theta(j))*[0:element_num-1]');p(j)=w'*a;endfigure;plot(theta,abs(p)),grid onxlabel('theta/radian')ylabel('amplitude')title('8阵元均匀线阵方向图')当来波方向为45度时，仿真图如下：8阵元均匀线阵方向图如下，来波方向为0度，20log（dB）随着阵元数的增加，波束宽度变窄，分辨力提高：仿真图如下：2.波束宽度与波达方向及阵元数的关系clcclear allclose allima=sqrt(-1);element_num1=16; %阵元数element_num2=128;element_num3=1024;lamda=0.03; %波长为0.03米d=1/2*lamda; %阵元间距与波长的关系theta=0:0.5:90;for j=1:length(theta);fai(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lamda/(element_num1*d));psi(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lamda/(element_num2*d));beta(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lamda/(element_num3*d)); endfigure;plot(theta,fai,'r',theta,psi,'b',theta,beta,'g'),grid onxlabel('theta');ylabel('Width in radians')title('波束宽度与波达方向及阵元数的关系')仿真图如下：3. 当阵元间距时，会出现栅瓣，导致空间模糊。

Capon波束形成器的仿真

Capon 波束形成器的仿真一、原理1、波束形成的定义波束形成就是从传感器阵列重构源信号，然后通过增加期望信源的贡献和抑制掉干扰源来实现的。

实际上波束形成可以看成空域（或时域）滤波器。

2、波束形成器的最佳权向量的推导首先，假定一天线阵列中个阵元的接收信号向量为()n x ，权向量为12[,,,]M w w w =w ，由于阵列的输出是对阵元的接收信号向量在个阵元上个分量的加权和，所以阵列的输出可写作：*1()()()MHmm m n n w x n ===∑y w x再令空间远场期望信号为()d t ，J 个干扰信号为()j i t ，1,2,,j J =，每个阵元上的加性白噪声为()k n t ，他们都具有相同的方差2n σ。

在这些假设条件下，第k 个阵元上的接收信号可以表示为1()()()()()()j Jk k d k i j k j x t a d t a i t n t θθ==++∑111212()()()()()()(),(),()()()()J d i i J M M x t d t n t x t i t n t a a a i t n t x t θθθ⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎡⎤=+⎣⎦⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦ 矩阵形式如下：1()()()()()()()()j Jd i j j t t t a d t a i t n t θθ==+=++∑x As n由以上可知：N 个快拍的波束形成器的输出为()()(1,2,,)H y t w x t t N ==输出的平均功率为2211222222111111()()()111()()()()()j N N H t t N J NNH Hd j i t j t t P w y t w x t N N w a d t i t w a w n t N N N θθ========⎡⎤=++⎢⎥⎣⎦∑∑∑∑∑∑当N →∞时，上式变为{}{}{}{}22222221()()()()()()()()jHHHJHHdj i nj P w E y t w E x t x t w wRwE d t w a E i t w a wθθσ=====++∑为了保证来自方向d θ的期望信号的正确接收，并完全抑制其它J 个干扰，关于权向量的约束条件应为：()1H d w a θ=和()0j H i w a θ=加上约束条件后{}222()()n P w E d t w σ=+波束形成器最佳权向量的确定以叙述为：在置零约束条件的约束下，求满足{}{}2min ()min HwwE y t w Rw =的权向量w 。

波束形成matlab程序

psi(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lamda/(element_num2*d));
beta(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lamda/(element_num3*d));
end
a=exp(imag*2*pi*d_lamda*sin(theta(j))*[0:element_num-1]'); p(j)=w'*a; end figure; plot(theta,abs(p)),grid on xlabel('theta/radian') ylabel('amplitude') title('8 阵元均匀线阵方向图 ')
当来波方向为 45 度时，仿真图如下：
8 阵元均匀线阵方向图如下，来波方向为 0 度， 20log（ dB ）
随着阵元数的增加，波束宽度变窄，分辨力提高：仿真图如下：
2. 波束宽度与波达方向及阵元数的关系
clc
clear all
close all
ima=sqrt(-1);
element_num1=16;
1.均匀线阵方向图
%8 阵元均匀线阵方向图，来波方向为 0 度 clc; clear all; close all; imag=sqrt(-1); element_num=8;% 阵元数为 8 d_lamda=1/2;% 阵元间距 d 与波长 lamda 的关系 theta=linspace(-pi/2,pi/2,200); theta0=0;% 来波方向 w=exp(imag*2*pi*d_lamda*sin(theta0)*[0:element_num-1]'); for j=1:length(theta)

1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。

Capon波束形成
阵列N=16,信号
0-30
θ︒
=，干扰为
160
θ︒=，
219
θ︒
=，
345
θ︒
=，干扰功率分别为：40dB,35dB,50dB。

Capon波束形成后的方向图和功率谱如下：
为了比较接收数据直接估计噪声协方差矩阵和利用干扰+噪声估计协方差矩阵的Capon 波束形成的差异，进行如下仿真：
可以看出利用干扰+噪声估计协方差矩阵的方向图性能较优于接收数据直接估计噪声协方差矩阵的方向图。

代码：
clc;
clear all ;
close all;
ima=sqrt(-1);
element_num=8; %阵元数
d_lamda=1/2; %阵元间距与波长的关系
theta=-90:0.5:90; %范围
theta0=-30; %来波方向
theta1=60; %干扰方向1
theta2=19; %干扰方向2
theta3=45; %干扰方向3
L=1000; %采样单元数
for i=1:L;
amp0=10*randn(1);%信号的幅度随机产生，保证信号之间是不相关的amp1=100*randn(1);%输入阵列的噪声
amp2=sqrt(10^3.5)*randn(1);%输入阵列的噪声
amp3=sqrt(10^5)*randn(1);%输入阵列的噪声
ampn=3;%噪声
x(:,i)=amp0*exp(ima*2*pi*1/2*sin(theta0*pi/180)*[0:element_num-1]')+...
amp1*exp(ima*2*pi*1/2*sin(theta1*pi/180)*[0:element_num-1]')+...
amp2*exp(ima*2*pi*1/2*sin(theta2*pi/180)*[0:element_num-1]')+...
amp3*exp(ima*2*pi*1/2*sin(theta3*pi/180)*[0:element_num-1]')+...
ampn*(randn(element_num,1)+ima*randn(element_num,1));
end
Rx=1/L*x* x';
R=inv(Rx);
steer=exp(ima*2*pi*1/2*sin(theta0*pi/180)*[0:element_num-1]');
w=R*steer/(steer'*R*steer);%Capon最优权矢量
for j=1:length(theta);
a=exp(ima*2*pi*d_lamda*sin(theta(j)*pi/180)*[0:element_num-1]');
f(j)=w'*a;
p(j)=1/(a'*R*a);
end
F=20*log10(abs(f)/(max(abs(f))));
P=20*log10(abs(p)/(max(abs(p))));%此处是功率的对数形式
figure;
% subplot(121)
plot(theta,F),grid on,hold on
plot(theta0,-80:0,'.')
plot(theta1,-80:0,'.')
plot(theta2,-80:0,'.')
plot(theta3,-80:0,'.')
xlabel('theta/o');
ylabel('F/dB');
title('Capon beamforming方向图') % axis([-90 90 -50 0]);
% subplot(122)
figure;
plot(theta,P),grid on
xlabel('theta/o');
ylabel('功率/dB');
title('Capon beamforming功率谱')。