作者：桂。
時間：2017-08-22 10:56:45
鏈接：http://www.cnblogs.com/xingshansi/p/7410846.html
前言
本文主要記錄常見的波束形成問題，可以說空間譜估計是波束形成基礎上發展而來，在系統論述空間譜之前，有必要分析一些Beamforming的基本特性。
一、波束形成模型
以均勻線陣為例：
按窄帶模型分析：
可以寫成矩陣形式：
其中
為方向矢量或導向矢量（Steering Vector），波束形成主要是針對各個接收信號X進行權重相加。
二、波束形成基本理論
A-波束形成
權重相加：
不同的波束形成，就是不同的權重W。
B-瑞利限
以均勻直角窗為例：
得出方向圖：
主瓣寬度正比於孔徑寬度的倒數：
因為孔徑的限制，造成波束寬度存在限制（不會無限制小），近而落在主瓣波束內部的兩個信號便會混在一起而分不清，這就存在瑞利限的問題。
直角窗主瓣寬度為：
其中λ為入射波長，theta1為入射角，Md為陣列孔徑。
C-常見窗函數
對於空間不同的陣列信號，類似采樣分析（空域采樣），自然可以加窗進行處理，不加窗可以認為是直角窗，另外也可以選擇漢明窗、hanning窗等等。
加窗可以改變波束寬度以及主瓣、副瓣等特性，可以借助MATLAB 的wvtool觀察不同窗函數特性。
N = 192; w = window(@blackmanharris,N); wvtool(w)
D-DFT實現
陣列的采樣間隔是相位信息：
這就類似於頻域變換，只不過這里的相位信息：對應的不是頻率，而是不同位置，可以看作空域的變換。
分別對陣列信號進行直接加權、加窗、DFT實現：
function x = StatSigGenerate(M,N,DOA,SNR,SignalMode,lambda,d)
ld = length(DOA);
if strcmp(SignalMode,'Independent')
st = randn(ld,N)+1j*randn(ld,N);
elseif strcmp(SignalMode,'Coherent')
st = [];
st1 = randn(1,N)+1j*randn(1,N);
for k = 1:ld
st = [st;st1];
end
end
st = st/sqrt(trace(st*st'/N)/ld);
nt = randn(M,N)+1j*randn(M,N);
nt = nt/sqrt(trace(nt*nt'/N)/M);
SNR = ones(1,ld)*SNR;
Amp = diag(10.^(SNR/20));
A = exp(1j*2*pi*[0:M-1]'*sind(DOA)*d/lambda);
x = A*Amp*st+nt;
end
主程序：
clc;clear all;close all
M = 32;
DOA = [-30 30];
SNR = 10;
theta = -90:.1:90;
len = length(theta);
SignalMode = 'Independent';
fc = 1e9;
c = 3e8;
lambda = c/fc;
d = lambda/2;
N = 100;%snap points
x = StatSigGenerate(M,N,DOA,SNR,SignalMode,lambda,d);
R_hat = 1/N*x*x';
output = zeros(3,len);
for i = 1:len
a = exp(1j*2*pi*[0:M-1]'*sind(theta(i))*d/lambda);
W = (inv(R_hat)*a)*(1./(a'*inv(R_hat)*a));
output(1,i) = mean(abs(W'*x),2);
output(2,i) = 1./(a'*inv(R_hat)*a);
output(3,i) = a'*x*ones(N,1);
end
output = abs(output);
output = output - repmat(mean(output.')',1,size(output,2));
output = output./repmat(max(output.')',1,size(output,2));
%plot
plot(theta,output(1,:),'k',theta,output(2,:),'r--',theta,output(3,:),'b');
legend('MVDR 波束','MVDR 譜','固定權重 波束');
對應結果圖：
E-自適應波束形成
直接相加也好、加窗也好，都是固定的權重系數，沒有考慮到信號本身的特性，所以如果結合信號本身去考慮就形成了一系列算法：自適應波束形成。
這類步驟通常是：
1）給定准則函數；
2）對准則函數進行求解。
准則常用的有：信噪比（snr）最大准則、均方誤差最小准則（MSE）、線性約束最小方差准則（LCMV）、最大似然准則（ML）等等；
求解的思路大體分兩類：1）直接求解，例如MVDR中的求解；2）也可以利用梯度下降的思想，如隨機梯度下降、批量梯度下降、Newton-raphson等方法，不再詳細說明。
以MVDR舉例：
這里采用直接求解的思路：
將求解的W帶入
即可得到波束形成。
F-柵瓣現象
柵瓣是一類現象，對應干涉儀就是相位模糊（相位超過2*pi），對應到Beamforming就是柵瓣問題，具體不再論述，給出現象（同樣的波束，在不同的位置分別出現）：
G-波束形成與空間譜
之前分析過MVDR的方法，得到的輸出（含有約束的最小均方誤差准則）為：
有時候也稱這個輸出為空間譜，其實就是|y2(t)|，但這個與MUSIC等算法的譜還不是一回事，只是有時候也被稱作空間譜，所以這里多啰嗦幾句，分析這個說法的來源。
已知N個采樣點的信號
，對其進行傅里葉變換：
進一步得到功率譜密度：
根據上文的分析：y(t)其實對應的就是空域變換（可借助DFT實現），類比於時頻處理中的頻域變換。而這里又可以看到頻域變換的平方/長度，對應就是功率譜，這是頻域的分析。
對應到空域，自然就是|y2(t)|/長度，對應空間譜，長度只影響比例關系，所以MVDR的最小方差輸出被稱作：空間譜也是合適的。
給出一個測試（這里如果），對比MVDR的y(t)、MVDR功率譜以及普通Beamforming的結果：
clc;clear all;close all
M = 32;
DOA = [-30 30];
SNR = 10;
theta = -90:.1:90;
len = length(theta);
SignalMode = 'Independent';
fc = 1e9;
c = 3e8;
lambda = c/fc;
d = lambda/2;
N = 100;%snap points
x = StatSigGenerate(M,N,DOA,SNR,SignalMode,lambda,d);
R_hat = 1/N*x*x';
output = zeros(3,len);
for i = 1:len
a = exp(1j*2*pi*[0:M-1]'*sind(theta(i))*d/lambda);
W = (inv(R_hat)*a)*(1./(a'*inv(R_hat)*a));
output(1,i) = mean(abs(W'*x),2);
output(2,i) = 1./(a'*inv(R_hat)*a);
output(3,i) = a'*x*ones(N,1);
end
output = abs(output);
output = output - repmat(mean(output.')',1,size(output,2));
output = output./repmat(max(output.')',1,size(output,2));
%plot
plot(theta,output(1,:),'k',theta,output(2,:),'r--',theta,output(3,:),'b');
legend('MVDR 波束','MVDR 譜','固定權重 波束');
對應結果：
如果將d = lambda/2;改為d = lambda/0.5;，自然就有了柵瓣：
