BYR Achieve · 镜像论坛

自己找了个程序。。估计有毛病。而自己是个菜鸟找不出。。请求大虾帮忙。。程序只有信道估计的图形。。而联合检测的图形没有。。误码率曲线也没有。。如何写程序使之出来呢。。还有就是程序的误码率太大了。如何修改程序呢。谢谢。。。望高手出现。。。附件(13.1KB) 源程序已附上。。。其中子程序加了%。。自己仿真时注意这个。。谢谢。。。。求高手帮忙。。。 % 程序File_C9:TD_JD.m clear; clc; %*********************************** %%仿真过程用到的变量 bit_len=88; %每个突发帧结构所包含的用户数据比特数 %TD-SCDMA每个突发帧的数据比特为704，扩频因子为16，采用QPSK调制，故每户数据比特数为704/16*2=88 N=44; %经过QPSK调制知道每个突发帧所包含的符号数 K=8; %用户个数 P=128; %基本midamble码的长度 Q=16; %扩频因子 W=P/K; %信道窗长度，为满足信道估计中的矩阵G是Toeplitz矩阵 %****************************************************************************** %%产生K个用户的随机数据 b=randint(K,bit_len); %随机产生K个用户的突发帧，实际系统中由信源产生 %******************************************************************************** %%进行QPSK调制 d=zeros(K,bit_len/2); for user=1:K for jj=1:bit_len/2 if b(user,2*jj-1)==0&&b(user,2*jj)==0 d(user,jj)=j; elseif b(user,2*jj-1)==0&&b(user,2*jj)==1 d(user,jj)=1; elseif b(user,2*jj-1)==1&&b(user,2*jj)==0 d(user,jj)=-1; elseif b(user,2*jj-1)==1&&b(user,2*jj)==1 d(user,jj)=-j; end end end %******************************************************************************* %%扩频 ovsf=[1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1; 1 1 1 1 1 1 1 1 -1 -1 -1 -1 -1 -1 -1 -1; 1 1 1 1 -1 -1 -1 -1 1 1 1 1 -1 -1 -1 -1; 1 1 1 1 -1 -1 -1 -1 -1 -1 -1 -1 1 1 1 1; 1 1 -1 -1 1 1 -1 -1 1 1 -1 -1 1 1 -1 -1; 1 1 -1 -1 1 1 1 1 -1 -1 1 1 -1 -1 1 1; 1 1 -1 -1 -1 -1 1 1 1 1 -1 -1 -1 -1 1 1; 1 1 -1 -1 -1 -1 1 1 -1 -1 1 1 1 1 -1 -1; 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1; 1 -1 1 -1 1 -1 1 -1 -1 1 -1 1 -1 1 -1 1; 1 -1 1 -1 -1 1 -1 1 1 -1 1 -1 -1 1 -1 1; 1 -1 1 -1 -1 1 -1 1 -1 1 -1 1 1 -1 1 -1; 1 -1 -1 1 1 -1 -1 1 1 -1 -1 1 1 -1 -1 1; 1 -1 -1 1 1 -1 -1 1 -1 1 1 -1 -1 1 1 -1; 1 -1 -1 1 -1 1 1 -1 1 -1 -1 1 -1 1 1 -1; 1 -1 -1 1 -1 1 1 -1 -1 1 1 -1 1 -1 -1 1]; %扩频码OVSF，Q=16时 spread=zeros(K,N*Q); %用扩频码扩频 for user=1:K spread(user,:)=kron(d(user,:),ovsf(user,:)); end %********************************************************************************* %%加扰 scram=[-1 1 -1 -1 -1 1 -1 -1 1 -1 1 1 -1 1 -1 -1]; w=ones(1,N); scram=kron(w,scram); s=zeros(size(spread)); for user=1:K s(user,:)=spread(user,:).*scram; end %******************************************************************************* %%形成突发结构 %突发帧结构由两个长度分别为352Chips的数据块， %一个长为144Chips的midamble和一个长为16Chips的保护间隔组成 %发送数据块E，对每个时隙，数据分为两个部分，中间放midamble码，最后放保护间隔 for user=1:K %总用户数K [mid_amble,GP]=midamble(K,P,user); %用户user的midamble码及保护间隔 len1=length(s(user,:))/2; %用户user的数据长度 len2=len1; len3=length(mid_amble); %用户user的midamble码长度 len=len1+len2+len3+16; %用户user的发送数据块长度 e=zeros(1,len); %用户user的数据块e for ii=1:len1 e(ii)=s(user,ii); %第一部分数据 e(ii+len1+len3)=s(user,ii+len1); %第二部分数据 end for ii=1:len3 e(ii+len1)=mid_amble(ii); %插入midamble码 end for ii=1:16 e(ii+len1+len3+len2)=GP(ii); end E(user,:)=e; %得到最后发送的包含midamble码的突发数据块E end %****************************************************************************************** %%debug D_E=sum(E(:,:)); R_DE=real(D_E); I_DE=imag(D_E); %**************************************************************** %%虚实分离 Re=real(E); Im=imag(E); %******************************************************************* %%脉冲成型 %设计升余弦滤波器 Delay=3; R=0.22; Fd=1; Fs=4; %生成升余弦滤波器 [ys,ts]=rcosine(Fd,Fs,'fir/sqrt',R,Delay); %Fd为输入序列的抽样率，Fs为升余弦滤波器的抽样率 %'fir/sqrt'表示滤波器为有限冲激的根升余弦 %R为滤波器参数Alpha,Delay为响应时延 %实部虚部分别通过滤波器 for user=1:K [y1(user,:),t1]=rcosflt(Re(user,:),Fd,Fs,'filter',ys); %length(y)=(length(x)+2*delay)*Fs/Fd [y2(user,:),t2]=rcosflt(Im(user,:),Fd,Fs,'filter',ys); send_data(user,:)=y1(user,:)+j*y2(user,:); end %******************************************************************** %%经过多径信道 %生成瑞利多径信道 Ts=0.02; %码片抽样时间间隔 Fmax=2; %最大多普勒频移 tau=[1,2,3,5]; %多径时延向量 pdb=[3,4,6,9]; %功率衰落向量 chan=rayleighchan(Ts,Fmax,tau,pdb); %必须满足Fmax<1/(10*Ts),chan为信道类 channel_hh=[1,zeros(1,243)]; %生成长度为244的冲激响应 channel_h=filter(chan,channel_hh); %将信道类转换成为信道响应向量 for jj=1:61 sample_h(jj)=channel_h(jj*4); %抽样得到长度为61的信道响应向量 end; for user=1:K channel_out(user,:)=conv(send_data(user,:),sample_h); %信号通过瑞利多径信道， %序列长度增加了w-1，w为信道窗长 end; channel_sum=sum(channel_out(:,:)); %****************************************************************************** %%经过高斯白噪声信道 %G_channel_sum=awgn(channel_sum,1,'measured'); %**************************************************************** %%接收端虚实分离 r_It=real(channel_sum); r_Qt=imag(channel_sum); %************************************************************************** %%匹配滤波 [yr1,tr1]=rcosflt(r_It,Fd,Fs,'filter/Fs',ys); %'filter/Fs'表示输入序列的抽样率等于滤波器抽样率为Fs [yr2,tr2]=rcosflt(r_Qt,Fd,Fs,'filter/Fs',ys); %length(y)=length(x)+(2*delay*Fs/Fd) ryr1=yr1.'; ryr2=yr2.'; %********************************************************************************** %%debug 抽样 for jj=0:N*Q+144+16+W-2 D_Isample(1,jj+1)=ryr1(1,25+jj*4); %抽样，从第25个开始每4个码片抽取1个，前24个码 %片为时延 D_Qsample(1,jj+1)=ryr2(1,25+jj*4); end D_Rresult=R_DE-D_Isample(1:864); D_Iresult=R_DE-D_Qsample(1:864); data=D_Isample+j*D_Qsample; %抽样得到的数据 %************************************************************************************************* %%信道估计 %Steiner估计器，该估计器只适合对FIR信道进行估计 %提取midamble码中的第W到Lm-1个 Lm=144; f_mid=N*Q/2+W; l_mid=N*Q/2+Lm-1; e=data(1,f_mid:l_mid); %构造G矩阵 basic=[1,-1,1,1,-1,-1,1,-1,1,-1,1,-1,1,1,-1,-1,-1,1,-1,-1,-1,-1,1,-1,-1,-1,-1,-1,1,1,1,1,... -1,1,1,1,1,1,-1,-1,1,-1,-1,-1,1,1,-1,1,1,1,1,-1,1,-1,1,1,1,1,1,1,1,-1,1,-1,... -1,1,1,-1,1,-1,-1,1,-1,1,-1,1,-1,-1,-1,-1,-1,1,-1,1,1,-1,-1,1,1,-1,-1,-1,-1,-1,-1,1,... 1,-1,1,1,1,1,-1,-1,1,1,-1,1,-1,-1,-1,-1,-1,-1,1,-1,1,-1,-1,-1,1,1,-1,-1,-1,-1,1,1]; %基本midamble码 %Lm is a constant,here is 144 L=Lm+(K-1)*W; %spread midamble length m=zeros(1,L); for ii=1:L if ii<=P m(1,ii)=(j^ii)*basic(1,ii); else m(1,ii)=m(1,ii-P); end end %产生G矩阵 G=zeros(P,P); %矩阵G是Toeplitz矩阵，也是一个循环矩阵 %其各列或各行由基本midamble循环移位得到 Position=0; for ii=1:P for jj=1:P Position=mod(P+ii-jj-1,P)+1; G(ii,jj)=m(Position); end end inv_G=inv(G); %求逆 hh=inv_G*e.'; %求得信道向量hh=[h1,h2,...hk] hh_s=hh; for user=1:K hhk(user,:)=hh_s((user-1)*W+1:user*W); %分别求出每个用户的信道向量 subplot(2,4,user); ph=1:1:W; figure(1) plot(ph,abs(hhk(user,:))); %画出每个用户的信道冲激响应，如果不加高斯白噪声，各个用户的信道冲激响应 %xlabel('信道冲激响应'）; %是相等的 %ylabel('响应窗长'）; end; % %**************************************************************************************************** %%联合检测 %求用户userde组合信道响应H及系统矩阵A A=zeros(N*Q+W-1,K*N); I_N=eye(N); for user=1:K scram=[-1 1 -1 -1 -1 1 -1 -1 1 -1 1 1 -1 1 -1 -1]; C_k=kron(I_N,(ovsf(user,:).*scram).'); %建立(N*Q,N)的块对角扩频矩阵 width=N*Q+W-1; H_k_widen=zeros(width,width); for ii=1:width indx_hh=1; col=ii; while (col>0)&&(indx_hh<=W) H_k_widen(ii,col)=hhk(user,indx_hh); col=col-1; indx_hh=indx_hh+1; end end H_k=H_k_widen(:,1:N*Q); %信道矩阵H_k为(N*Q+W-1,N*Q)的矩阵 B_k=H_k*C_k; %矩阵B_k为(N*Q+W-1,N)的矩阵 A(:,(user-1)*N+1:user*N)=B_k; %矩阵A=[B_k,B_K-1,...B_1]为(N*Q+W-1,K*N)的矩阵 end % figure(2) % plot(ii,user); %分离出接收的数据 x1=data(1:N*Q/2); %第一数据块 x2=data(N*Q/2+145:N*Q+144+15); %第二数据块加信道窗响应 xx=[x1,x2]; %根据迫零算法，接收端的输出信号估计为 a_A=A'; AA=inv(a_A*A); d_estimate=AA*a_A; dd_estimate=d_estimate*(xx.'); d_restore=dd_estimate.'; %***************************************************************************************************************** %%QPSK逆映射 for user=1:K dd_restore(user,:)=d_restore((user-1)*N+1:user*N); end; bb=zeros(K,88); for user=1:K for jj=1:N for kk=1:4 test(kk)=abs(dd_restore(user,jj)-j^kk); %求出接收信号与QPSK判决信号的欧氏距向量 end; test_result=min(test); %加入最小欧氏距判决，并进行QPSK逆映射 if test_result==test(1) bb(user,2*jj-1)=0; bb(user,2*jj)=0; elseif test_result==test(2) bb(user,2*jj-1)=1; bb(user,2*jj)=0; elseif test_result==test(3) bb(user,2*jj-1)=1; bb(user,2*jj)=1; elseif test_result==test(4) bb(user,2*jj-1)=0; bb(user,2*jj)=1; end; end; end; % figure(2) % plot(jj,kk); %************************************************************************************************** %%求误码率 error=bb-b; for user=1:K count(user)=0; for jj=1:bit_len if error(user,jj)~=0 count(user)=count(user)+1; end; end; end; count figure(2) plot(user,jj); %********************************************************************************** % 程序结束 % % 子程序： % % function[mid_amble,GP]=midamble(K,P,user) % ************************************************************************************ % 输入参数： % K：用户数 % P：基本midamble码长度 % user:使用所产生midamble码的用户 % 输出参数： % mid_amble:所产生的midamble码 % GP：所产生的保护间隔 % ************************************************************************** % % 生成复数midamble码 % W=P/K; %W为冲激窗长 % Lm=144; %Lm为生成的midamble码长 % basic=[1,-1,1,1,-1,-1,1,-1,1,-1,1,-1,1,1,-1,-1,-1,1,-1,-1,-1,-1,1,-1,-1,-1,-1,-1,1,1,1,1,... % -1,1,1,1,1,1,-1,-1,1,-1,-1,-1,1,1,-1,1,1,1,1,-1,1,-1,1,1,1,1,1,1,1,-1,1,-1,... % -1,1,1,-1,1,-1,-1,1,-1,1,-1,1,-1,-1,-1,-1,-1,1,-1,1,1,-1,-1,1,1,-1,-1,-1,-1,-1,-1,1,... % 1,-1,1,1,1,1,-1,-1,1,1,-1,1,-1,-1,-1,-1,-1,-1,1,-1,1,-1,-1,-1,1,1,-1,-1,-1,-1,1,1,]; % 基本midamble码，长度为P % L=Lm+(K-1)*W; %进行周期扩展 % m=zeros(1,L); % for ii=1:L %按照“标准”生成复数midamble码 % if ii<=P % m(1,ii)=(j^ii)*basic(1,ii); % else % m(1,ii)=m(1,ii-P); % end % end % 生成用户所使用的midamble码 % % mid_amble=zeros(1,Lm); % for ii=1:Lm % mid_amble(1,ii)=m(1,ii+(K-user)*W); % end % 产生保护间隔Gp % Gp=2*random('bino',1,0.5,1,16)-1; % for k=1:16 % if(mod(k,2)~=0) % GP(k)=GP(k)*j; % end % end % 程序结束

matalb信道估计及联合检测。。急救！！！！