新建 文本文档.txt

信号分选中关于信号的盲分离

clc;clear all;close all; 
 
%%%%%%%%%%%%%%  读入原始图像，混合，并输出混合图像  %%%%%%%%%%%%%%%%%% 
 
% 读入混合前的原始图片并显示 
t=0:1/100:9; 
I1=sin(t);  
I2=randn(1,901); 
I3=square(4*t); 
subplot(4,3,1),plot(I1),title('输入信号1'), 
subplot(4,3,2),plot(I2),title('输入信号2'), 
subplot(4,3,3),plot(I3),title('输入信号3'), 
 
% 将其组成矩阵 
S=[I1;I2;I3];                          % 图片个数即为变量数，图片的像素数即为采样数 
                                       % 因此S_all是一个变量个数＊采样个数的矩阵 
 
Sweight=randn(size(S,1));              % 取一随机矩阵，作为信号混合的权矩阵 
MixedS=Sweight*S;                      % 得到三个混合信号矩阵 
 
% 将混合矩阵重新排列并输出 
subplot(4,3,4),plot(MixedS(1,:)),title('混合信号1'), 
subplot(4,3,5),plot(MixedS(2,:)),title('混合信号2'), 
subplot(4,3,6),plot(MixedS(3,:)),title('混合信号3'), 
 
MixedS_bak=MixedS;                      % 将混合后的数据备份，以便在恢复时直接调用 
%%%%%%%%%%%%%%%%%%%%%%%%%%  标准化  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
MixedS_mean=zeros(3,1); 
for i=1:3 
    MixedS_mean(i)=mean(MixedS(i,:)); 
end                                        % 计算MixedS的均值 
 
for i=1:3 
    for j=1:size(MixedS,2) 
        MixedS(i,j)=MixedS(i,j)-MixedS_mean(i); 
    end 
end                                        
 
%%%%%%%%%%%%%%%%%%%%%%%%%%%  白化  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
 
MixedS_cov=cov(MixedS');                    % cov为求协方差的函数 
[E,D]=eig(MixedS_cov);                      % 对图片矩阵的协方差函数进行特征值分解 
Q=inv(sqrt(D))*(E)';                        % Q为白化矩阵 
MixedS_white=Q*MixedS;                      % MixedS_white为白化后的图片矩阵 
IsI=cov(MixedS_white');                     % IsI应为单位阵             
 
%%%%%%%%%%%%%%%%%%%%%%%%　FASTICA算法  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
X=MixedS_white;                            % 以下算法将对X进行操作 
[VariableNum,SampleNum]=size(X); 
numofIC=VariableNum;                       % 在此应用中，独立元个数等于变量个数 
B=zeros(numofIC,VariableNum);              % 初始化列向量w的寄存矩阵,B=[b1  b2  ...   bd] 
for r=1:numofIC                            % 迭代求取每一个独立元 
    i=1;maxIterationsNum=100;              % 设置最大迭代次数（即对于每个独立分量而言迭代均不超过此次数） 
    IterationsNum=0; 
    b=rand(numofIC,1)-.5;                  % 随机设置b初值 
    b=b/norm(b);                           % 对b标准化 
    while i<=maxIterationsNum+1 
        if i == maxIterationsNum           % 循环结束处理 
            fprintf('\n第%d分量在%d次迭代内并不收敛。', r,maxIterationsNum); 
            break; 
        end 
        bOld=b;                       
        a2=1; 
        u=1; 
        t=X'*b; 
        g=(exp(2.*t)-1)./(exp(2.*t)+1); 
        dg=4*exp(2.*t)./(exp(2.*t)+1).^2; 
        b=((1-u)*t'*g*b+u*X*g)/SampleNum-mean(dg)*b; 
                                            % 核心公式，参见理论部分公式2.52 
        b=b-B*B'*b;                         % 对b正交化 
        b=b/norm(b);  
        if abs(abs(b'*bOld)-1)<1e-9         % 如果收敛，则保存b 
               B(:,r)=b; 
             break; 
         end 
        i=i+1;         
    end 
end 
 
%%%%%%%%%%%%%%%%%%%%%%%%%  数据复原并构图  %%%%%%%%%%%%%%%%%%%%%%%% 
 
ICAedS=B'*Q*MixedS_bak;                     % 参见理论部分公式2.55 
 
% 将混合矩阵重新排列并输出 
subplot(4,3,7),plot(ICAedS(1,:)),title('ICA输出信号1'), 
subplot(4,3,8),plot(ICAedS(2,:)),title('ICA输出信号2'), 
subplot(4,3,9),plot(ICAedS(3,:)),title('ICA输出信号3'),