css.m

宽带信号测向的几种直接处理算法

function P=CSS(data,angles,fl,fh,fs,D,theta0,x);
% 宽带聚焦高分辨方法——CSS。
%     P=CSS(data,angles,fl,fh,fs,D,theta0,x);
% 
%输入参数：
%   data    -- 阵列接收的数据，每行为一次快拍，每列为一个阵元输出；
%   angles   -- 方位搜索范围向量,如[0:30];
%   fh      -- 上边频；
%   fl      -- 下边频；
%   fs      -- 采样频率；
%   D       -- 信号子空间有效维数；
%   theta0  -- 聚焦方位;
%   x       -- 线列阵各阵元的坐标(米)，如果给一个整数，表示x元的均匀线列阵阵元间距为中心频率的半波长,
%              如 [0:0.05:0.3]。
% 
% 输出参数：
%     P    -- 空间角谱。
%
%Example:
%  PCSS=CSS(data,angles,fl,fh,fs,dcss,theta0,x);
%
%See also SOM, BASS, FBASS.

% Reference:
%   Wang and M. Kaveh, "Coherent siganl-subspace processing for the detectino and
%     estimation of angles of arrival of multiple wide-band sources",
%     IEEE Trans. ASSP, vol.33, pp.823-831, 1985.
%

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    % Copyright(c) 2006-2026 by Zhang Lijie %
    % Directed by Zhang Qunfei              %
    % $Version 2.0$, $Date: 2007/3/30$      %
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% Built-out function

c=1490;             %声速
%%确保x为阵元坐标或阵元数目
if(length(x)==1)
    if x>1 && x<=100 && x==fix(x)
        f0=(fl+fh)/2;  lamda=c/f0;
        x=(0:x-1)*lamda/2;
    else
        error('阵元数x必须为2～100的整数！')
    end
end
%%对各列进行fft变换
[N,M]=size(data);  % N=10240;  %FFT点数超过数据点数没有意义
if nargin<5 || M>N,  help(mfilename),  end
data=fft(data,N);
data=data(1:fix(N/2),:);
%取出'fl'kHz－'fh'kHz的谱线
nl=fix(fl*N/fs);
nh=fix(fh*N/fs);
data=data(nl:nh,:);  %每行是特定频率的快拍，每列是一个阵元的频谱
f=[nl*fs/N:fs/N:nh*fs/N];
%%Broad-Band Covariance Matrix Estimation
[J,M]=size(data);
f0=(fl+fh)/2;
a0=exp(-1*j*2*pi*f0/c*x*sin(pi/180*theta0)).';
R=zeros(M);
for i=1:J
    ai=exp(-1*j*2*pi*f(i)/c*x*sin(pi/180*theta0)).';
    T=diag(a0./ai);
    X=data(i,:).';
    rr=X*X';
    R=R+T*rr*T';
end
RR=R/J^2;
%%空间谱估计
[V d]=eig(RR);
[U I]=sort(abs(diag(d)));
V=V(:,I);
En=V(:,1:M-D);
P=[];
for ang=angles
    A=exp(-1*j*2*pi*f0/c*x*sin(pi/180*ang)).';
    p=1/abs(A'*En*En'*A);
    P=[P,p];
end