基于时序结构的盲信号分离

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Copyright (c) by JV Stone, University of Sheffield. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% clear all; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % SET UP %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% global demo_id; %声明全局变量% % Set rand number seeds. seed=9;randn('state',seed);rand('state',seed); % randn ：产生标准正态分布的随机数或矩阵的函数， % RANDN('state',J) J 为整数，设置发生器到J阶状态 demo_id = 2; % Used by get_data, set to 1 or 2, default=2.. num_sources = 3; num_mixtures = num_sources; num_samples = 5000; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % GET DATA. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % sources = (num_samples x num_sources), % mixtures = (num_samples x num_mixtures), % A = (num_sources x num_sources) mixing matrix, mixtures=sources*A; % Get zero mean, unit variance sources. [sources mixtures A] = get_data(num_samples); % one mixture per column. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % COMPUTE V AND U. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Set short and long half-lives. shf = 1; lhf = 900000; max_mask_len= 500; n = 8; % n = num half lives to be used to make mask. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Get masks to be used to find (x_tilde-x) and (x_bar-x) % Set mask to have -1 as first element, % and remaining elements sum to unity. % Short-term mask. h=shf; t = n*h; lambda = 2^(-1/h); temp = [0:t-1]'; lambdas = ones(t,1)*lambda; mask = lambda.^temp; mask(1) = 0; mask = mask/sum(abs(mask)); mask(1) = -1; s_mask=mask; s_mask_len = length(s_mask); % Long-term mask. h=lhf;t = n*h; t = min(t,max_mask_len); t=max(t,1); lambda = 2^(-1/h); temp = [0:t-1]'; lambdas = ones(t,1)*lambda; mask = lambda.^temp; mask(1) = 0; mask = mask/sum(abs(mask)); mask(1) = -1; l_mask=mask; l_mask_len = length(l_mask); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Filter each column of mixtures array. S=filter(s_mask,1,mixtures); L=filter(l_mask,1,mixtures); % Find short-term and long-term covariance matrices. U=cov(S,1); V=cov(L,1); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Find eigenvectors W and eigenvalues d. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% [W d]=eig(V,U); W=real(W); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Recover source signals. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ys = mixtures*W; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % PLOT RESULTS %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Rescale ys to zero-mean and unit variance, for display purposes. temp=repmat(mean(ys),[num_samples 1]'); ys=ys-temp; temp=repmat(std(ys,1),[num_samples 1]'); ys=ys./temp; % Plot original sources ... figure(1);clf; nn=1000; if demo_id==1 imin=100; else imin=1; end; imin=1;imax=imin+nn-1; figure(1); for i=1:num_sources subplot(num_sources,1,i); plot(mixtures(imin:imax,i)/std(mixtures(imin:imax,i),1),'k.'); axis off; end; subplot(num_sources,1,1); title('Signal Mixtures'); % Plot corresponding recovered sources. a=[sources ys]; c=corrcoef(a); rs=c(1:num_sources,num_sources+1:num_sources*2); fprintf('Correlations between sources and recovered signals ...\n'); abs(rs) rs_bin=abs(rs)>0.9; figure(2); for i=1:num_sources subplot(num_sources,1,i); hold off; plot(-1+sources(imin:imax,i),'k'); hold on; % Find corresponding sources/ys via correlation matrix rs_bin. j=find(rs_bin(i,:)); s=sign(rs(i,j)); plot(2+s*ys(imin:imax,j),'r--'); hold off; axis off; end; subplot(num_sources,1,1); title('Source signal (bottom trace) and recovered signal (top trace)'); % Play sounds if demo_id==1 fprintf('Playing mixtures ...\n');pause(2); for i=1:num_sources soundsc(mixtures(:,i));pause(1); end; fprintf('Playing recovered signals ...\n');pause(2); for i=1:num_sources soundsc(ys(:,i));pause(1); end; end; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%