clear?all;
close?all;
clc;
%%?signal?model
source_number=16;%信元數
sensor_number=10;%陣元數
N_x=10240;?%信號長度
snapshot_number=N_x;%快拍數
l=0.15;%信號波長??
d=l/2;%陣元間距
snr=20;%信噪比
fs=pi/6;
%w1=unifrnd（121source_number）;
w1=1:1/（source_number-1）:2;
w=w1*fs;%信號頻率
w=w.‘;
S=sqrt（10.^（snr/10））*exp（j*w*[0:N_x-1]）;%仿真信號non-Gussian?signal
theta_range?=?50;
source_doa?=?（-theta_range:2*theta_range/（source_number-1）:theta_range）;
pos=[0:sensor_number-1]‘;
%pos=[1?2?3?6?9]‘;
A=exp（-j*pos*d*2*pi*sin（source_doa*pi/180）/l）;
x=A*S+（1/sqrt（2））*（randn（sensor_numberN_x）+j*randn（sensor_numberN_x））;%加了高斯白噪聲后的陣列接收信號
%%?fourth-cumulant
B=zeros（sensor_number^2source_number）;?
C_s=zeros（source_numbersource_number）;
for?a=1:source_number
????B（:a）=kron（A（:a）conj（A（:a）））;
end
for?b=1:source_number
????for?c=1:N_x
????????%?第一項
????????s_temp1=S（bc）*conj（S（bc））*conj（S（bc））*S（bc）;
????????%?第二項
????????s_temp2=S（bc）*conj（S（bc））;
????????s_temp3=conj（S（bc））*S（bc）;
????????%?第三項
????????s_temp4=S（bc）*S（bc）;
????????s_temp5=conj（S（bc））*conj（S（bc））;
????end
????s_temp1=s_temp1/N_x;
????s_temp2=s_temp2/N_x;
????s_temp3=s_temp3/N_x;
????s_temp4=s_temp4/N_x;
????s_temp5=s_temp5/N_x;
????????%?總項
????s_b=s_temp1-s_temp2*s_temp3-s_temp4*s_temp5-s_temp2*s_temp3;
????C_s（bb）=s_b;
end
C4=B*C_s*B‘;
%%?MUSIC
[VD]=eig（-C4）;%對C4進行特征分解
[lambda?index]=sort（diag（D））;
En=V（:index（1:（sensor_number^2-source_number）））;
disp（En）;
searching_doa=-90:0.5:90;%線陣的搜索范圍為-90~90度
?for?i=1:length（searching_doa）
???a_theta=exp（-j*（0:sensor_number-1）‘*2*pi*d*sin（pi*searching_doa（i）/180）/l）;
???b_theta=kron（a_thetaconj（a_theta））;
???Pmusic（i）=1./（norm（b_theta‘*En））.^2;
?end
?figure;
plot（searching_doa10*log（Pmusic））;
xlabel（‘入射角/degree‘）;
ylabel（‘功率譜/dB‘）;
title（‘四階累積量MUSIC估計‘）;
grid?;

?屬性????????????大小?????日期????時間???名稱
-----------?---------??----------?-----??----
?????文件????????1920??2020-03-09?21:26??四階累積量\fourth_Cum_MUSIC.m
?????文件??????500160??2020-03-09?21:36??四階累積量\基于四階累積量特征分解的空間譜估計測向方法_魏平.pdf
?????目錄???????????0??2020-03-09?21:45??四階累積量\