以下是我写的一个2发3收MIMO-FIR信道SE-CMA盲均衡算法的程序,用的是基于互信息的盲源分离方法,但是输出后还是锁定到一根发射天线,请高手指教一下,谢谢!!!
N=1000; % 信号序列长度
K=3; % 滤波器长度
snr=10;
delta=0.0002;
sigpower=1;
d=1; % 欧式空间距离
M=4; % QAM进制数
a=0.3;
mapping=[-d d;-d -d;d -d;d d];
% H[0]=[-1.9522 -0.5706;-0.5666 0.4246;-1.1293 0.7666]
% H[1]=[1.0691 -1.8841;-0.7926 0.0598;0.3569 -0.2744]
% 2发3收的MIMO-FIR信道
for i=1:N
temp1=rand;
dsource1(i)=1+floor(M*temp1);
temp2=rand;
dsource2(i)=1+floor(M*temp2);
end;
for i=1:N
qam_sig1(i,:)=mapping(dsource1(i),:);
qam_sig2(i,:)=mapping(dsource2(i),:);
end; % 信号映射
r1=qam_sig1;
rr1=r1(:,1)+j*r1(:,2);
r2=qam_sig2;
rr2=r2(:,1)+j*r2(:,2);
rr0(:,1)=rr1;
rr0(:,2)=rr2;
rr=rr0'; % 发送复数信号
for i=1:N
ya(1,i)=filter([-1.9522 1.0691],1,rr(1,i))+filter([-0.5706 -1.8841],1,rr(2,i));
ya(2,i)=filter([-0.5666 -0.7926],1,rr(1,i))+filter([0.4246 0.0598],1,rr(2,i));
ya(3,i)=filter([-1.1293 0.7666],1,rr(1,i))+filter([0.7666 -0.2744],1,rr(2,i));
end; % 序列通过2发3收的MIMO-FIR信道
yout(1,:)=awgn(ya(1,:),snr,sigpower);
yout(2,:)=awgn(ya(2,:),snr,sigpower);
yout(3,:)=awgn(ya(3,:),snr,sigpower);
ya1=yout(1,:);
ya2=yout(2,:);
ya3=yout(3,:);
g1=[0 0 0 1 1 1 0 0 0];
g2=[0 0 0 1 1 1 0 0 0]; % 均衡器抽头系数初值
for k=1:N-K+1 % 迭代次数,998次
y_k=[ya1(k+K-1) ya2(k+K-1) ya3(k+K-1) ya1(k+K-2) ya2(k+K-2) ya3(k+K-2) ya1(k) ya2(k) ya3(k)];
z_k1=g1*y_k';
z_k2=g2*y_k'; % 均衡器输出
yt1(k)=z_k1;
yt2(k)=z_k2;
H=[0 0 0 0 -1.9522 -0.5706 1.0691 -1.8841;0 0 0 0 -0.5666 0.4246 -0.7926 0.0598;0 0 0 0 -1.1293 0.7666 0.3569 -0.2744;0 0 -1.9522 -0.5706 1.0691 -1.8841 0 0;0 0 -0.5666 0.4246 -0.7926 0.0598 0 0;0 0 -1.1293 0.7666 0.3569 -0.2744 0 0;-1.9522 -0.5706 1.0691 -1.8841 0 0 0 0;-0.5666 0.4246 -0.7926 0.0598 0 0 0 0;-1.1293 0.7666 0.3569 -0.2744 0 0 0 0];
% MIMO-FIR信道卷积矩阵
% [0 0 H(0) H(1);0 H(0) H(1) 0;H(0) H(1) 0 0]
s(1,:)=g1*H;
s(2,:)=g2*H;
c11=[s(1,1) s(1,3) s(1,5) s(1,7)];
c12=[s(1,2) s(1,4) s(1,6) s(1,8)];
c21=[s(2,1) s(2,3) s(2,5) s(2,7)];
c22=[s(2,2) s(2,4) s(2,6) s(2,8)];
IT2(k)=(abs(s(1,2))^2+abs(s(2,2))^2+abs(s(1,4))^2+abs(s(2,4))^2+abs(s(1,6))^2+abs(s(2,6))^2+abs(s(1,8))^2+abs(s(2,8))^2-max([abs(s(1,2))^2 abs(s(2,2))^2 abs(s(1,4))^2 abs(s(2,4))^2 abs(s(1,6))^2 abs(s(2,6))^2 abs(s(1,8))^2 abs(s(2,8))^2]))/max([abs(s(1,2))^2 abs(s(2,2))^2 abs(s(1,4))^2 abs(s(2,4))^2 abs(s(1,6))^2 abs(s(2,6))^2 abs(s(1,8))^2 abs(s(2,8))^2]);
% 均衡后ISI和MUI的总和,以第2个均衡器输出为例
dj1=4*sign(abs(z_k1)-1)*sign(z_k1)*conj(y_k)+a*sign(z_k1*conj(z_k2))*conj(z_k2)*conj(y_k); % 代价函数的梯度
g1=g1-0.5*delta*dj1; % 抽头系数迭代更新
dj2=4*sign(abs(z_k2)-1)*sign(z_k2)*conj(y_k)+a*sign(z_k2*conj(z_k1))*conj(z_k1)*conj(y_k); % 代价函数的梯度
g2=g2-0.5*delta*dj2; % 抽头系数迭代更新
end;
figure(1);
plot(real(yt1),imag(yt1),'*');
figure(2);
plot(real(yt2),imag(yt2),'*');
figure(3);
subplot(2,2,1),stem(abs(c11));
subplot(2,2,2),stem(abs(c12));
subplot(2,2,3),stem(abs(c21));
subplot(2,2,4),stem(abs(c22));
figure(4);
plot(IT2);