K=8; % 信道非零抽头稀疏，即信道稀疏度；
N=128; % 信道长度
M=64; % 输出序列长度，即测量数(M>=K*log(N/K),至少40,但有出错的概率)
%1.产生2*2稀疏多径信道，共4条路径，每条路径稀疏度为K=8；
h=zeros(4,N);
p=zeros(1,K);
for i=1:4
for n=1:K
p(1,n)=randint(1,1,[1 N]);
h(i,p(1,n))=2*rand(1,1)-1;
end
end
%2.产生测量矩阵，两种方法：toeplitz矩阵与随机高斯矩阵；
%toeplitz矩阵
% X1=toeplitz_1(M,N);
% X2=toeplitz_1(M,N);
%随机高斯矩阵
X1=randn(M,N);
X2=randn(M,N);
%3.获得信道冲击响应观测值，暂时未加入噪声；
X=zeros(64,256);
X(:,1:128)=X1;
X(:,129:256)=X2;
y1=X1*h(1,:)'+X2*h(3,:)';
y2=X1*h(2,:)'+X2*h(4,:)';
%4.重构算法，预计使用传统的信道估计方法（LS,MMSE），及压缩感知估计（BP,OMP,CoSaMP）作图并对比
h1=omp(y1,X,2*N);
h11=h1(1,1:128);
h21=h1(1,129:256);
h2=omp(y2,X,2*N);
h12=h2(1,1:128);
h22=h2(1,129:256);
subplot(4,2,1);plot(h(1,:));title('发射天线1到接收天线1的信道冲激响应');
subplot(4,2,2);plot(h11);title('omp重构的响应');
subplot(4,2,3);plot(h(2,:));title('发射天线1到接收天线2的信道冲激响应');
subplot(4,2,4);plot(h12);title('omp重构的响应');
subplot(4,2,5);plot(h(3,:));title('发射天线2到接收天线1的信道冲激响应');
subplot(4,2,6);plot(h21);title('omp重构的响应');
subplot(4,2,7);plot(h(4,:));title('发射天线2到接收天线2的信道冲激响应');
subplot(4,2,8);plot(h22);title('omp重构的响应');
figure; %1发1收的原始与重构信号
hold on;
plot(h11,'k.-') % 重建信号
plot(h(1,:),'r*') % 原始信号
legend('Recovery','Original')
