📄 multi_awgn_vita_rsencoded_bis_all_wdiv2.m
字号:
%function pb=multi_awgn_vita(Eb_to_Nj_in_dB,Eb_to_No_in_dB,BPH,number_of_states,Q)
%VITERBI This program simulates the Viterbi sequnce decoding of the
% differential frequency hopping system with (15,7) RS
% pre-encoded and diversity of level 2
%
% Eb_to_Nj_in_dB is the signal-to-jamming ratio given in dB
%
% Eb_to_No_in_dB is the signal-to-noise ratio given in dB
%
% BPH is the number of bits transmitted by one hop
%
% number_of_states is the number of states in the DFH
% trellis,corresponding to the right L stages of the DFH encoding
% shift register
%
% Q is the number of jamming tones in the DFH bandwidth
%
% 这个程序是正确的,带(N,K,D)=(15,7,9)的 RS 纠错码,编码后在后面只补上L个0,即符号流总长度是2N+L而不是2(N+L),
%
%
%
%
%********系统参数********%
% Eb_to_Nj_in_dB=11;
Eb_to_Nj_in_dB=[0,1];
Eb_to_No_in_dB=13.35;
BPH=2;
number_of_states=16;
Q=2; % 干扰音的个数
N=700; % 每次符号流长度
times=5; % 重复做500次
PB_all=zeros(length(Eb_to_Nj_in_dB),1);
PS_all=zeros(length(Eb_to_Nj_in_dB),1);
vita_symbol_err=zeros(length(Eb_to_Nj_in_dB),times);% 维特比译码后, RS 译码之前2^BPH=4进制符号错误个数计数器
num_of_err=zeros(length(Eb_to_Nj_in_dB),times);% 最终 RS 码译码输出的2^M=16进制符号错误个数计数器
M=4; % 伽罗华域GF(2^m)中m=4,即GF(16)域
% 4位二进制序列与GF(16)域中1个16进制GF符号对应转换关系表
binary_table=[0 0 0 0;
0 0 0 1;
0 0 1 0;
0 0 1 1;
0 1 0 0;
0 1 0 1;
0 1 1 0;
0 1 1 1;
1 0 0 0;
1 0 0 1;
1 0 1 0;
1 0 1 1;
1 1 0 0;
1 1 0 1;
1 1 1 0;
1 1 1 1];
% 2位四进制序列与GF(16)域中1个16进制GF符号对应转换关系表
quarter_table=[0 0;
0 1;
0 2;
0 3;
1 0;
1 1;
1 2;
1 3;
2 0;
2 1;
2 2;
2 3;
3 0;
3 1;
3 2;
3 3];
decimal_table=[0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15]';
GF16_table=gf(decimal_table,M);% GF(16)域中GF符号表
%************************%
WTbarpp=waitbar(0,'Outer SNRs Loop:Please wait');
for pp=1:length(Eb_to_Nj_in_dB)
fanout=2^BPH; % DFH的扇出系数
Eb_to_Nj=10^(Eb_to_Nj_in_dB(pp)/10); % 比值形式的Eb/Nj
Eb_to_No=10^(Eb_to_No_in_dB/10); % 比值形式的Eb/No
L=floor(log(number_of_states)/log(fanout)); % 编码移位寄存器的长度为L+1,最右边L级是其状态位,与网格图中的状态一一对应(注意并不是与跳频频点一一对应)
% pb=zeros(1,times);
WTbar=waitbar(0,'SNR inside loop:please wait for a moment...');
for rep=1:times
% source=[randint(1,N,fanout),zeros(1,L)]; % 信息源:(注意不是10进制的,而是fanout进制的)随机符号流,最后补上L个0符号,使移位寄存器的状态清零
%***************信源部分(信源是长度为N的16进制符号流)*********************%
source_sequence=randint(1,N,2^M);
source_matrix=reshape(source_sequence,7,N/7); % 特别注意reshape函数是从矩阵A中按列的顺序取出元素,返回一个M*N的矩阵,但是注意是把元素按先填第1列,然后填第2列……的顺序填到新矩阵中去的，不是按行填写
source=source_matrix'; % reshape之后得到的是一个7行100列的矩阵,转置后得到的是一个100行7列的矩阵source,矩阵source的第一行存放的是source_sequence
% 序列的第1到7个符号,第二行存放的是source_sequence序列的第8到14个符号,依此类推
%*************************************************************************%
% %***********************信源部分************************%
%
% % 将信息符号流转化成二进制信息比特流
% dsource=zeros(1,N*BPH);
% if(BPH~=1)
% for i=1:N
% dsource((i-1)*BPH+1:i*BPH)=deci2change(source(i),BPH,2); %%%%%%%%%%%%其实这里是把 2^BPH 进制的符号source(i)转换成BPH个二进制比特,前面产生的符号流source=randint(1,N,fanout)
% %产生的符号流原来本身也就是fanout=2^BPH进制的,这正好相吻合,所以这里并不是把十进制的符号source(i)转化为
% %BPH位二进制比特,而是把2^BPH进制的符号source(i)(显然source=randint(1,N,fanout)表明source(i)
% %的取值范围被限定在0到fanout-1之间)转化为BPH位二进制比特,
% end
% else
% dsource=source(1:N);%*****
% end
%
% %******************************************************%
% % *******************差错控制编码部分(8,4)码**********************%
% % G=[1 0 0 0 1 0 1 1;
% % 0 1 0 0 1 1 1 0;
% % 0 0 1 0 1 1 0 1;
% % 0 0 0 1 0 1 1 1];
% G=[1 0 0 0 1 0 1 1;0 1 0 0 1 1 0 1;0 0 1 0 0 1 1 1;0 0 0 1 1 1 1 0]; % (8,4,4)扩展汉明码的生成矩阵G
% % H=[1 1 1 1 1 1 1 1;1 1 0 1 1 0 0 0;0 1 1 1 0 1 0 0;1 0 1 1 0 0 1 0]; % (8,4,4)扩展汉明码的一致校验矩阵H
% % trt=syndtable(H); % (8,4,4)扩展汉明码的伴随式译码表
% dsource_coded=encode(dsource,8,4,'linear',G)';
% % source_coded0=zeros(1,2*N);
% % *****************************************************************%
%****************差错控制编码部分(N,K,Dmin)=(15,7,9)的 RS 码*****************%
msg_GF=gf(source,M); %把信息符号转化为GF域的符号
source_coded_GF=rsenc(msg_GF,15,7); % RS编码后,source_coded_GF是一个100行15列的矩阵,矩阵每一行是一个 RS 码的码字
%****************************************************************************%
% % *********************交织部分(块交织)(基于比特的交织)****************%
% % 仅适用于BPH=2且N=1000时
% interleave_outcome=zeros(1,2*BPH*N);
% A1=zeros(200,200);
% for i=1:2*BPH*N
% A1(i)=dsource_coded(i);
% end
% A=A1.';
% for i=1:2*BPH*N
% interleave_outcome(i)=A(i);
% end
% %************************************************%
% % *******************随机交织************************%
% [interleave_outcome,alpha]=interleave(dsource_coded); % dsource_coded 长度为(n/k)*N*BPH=2*N*BPH
%
% % ***************************************************%
%**************交织部分(基于16进制GF(16)域符号(即RS码元)的交织)***************%
codedsequence=zeros(1,1500);
codedsequence_GF=gf(codedsequence,M);
for i=1:1500
codedsequence_GF(i)=source_coded_GF(i); % 把RS编码后的码字矩阵转化为一串码字序列,按列读出,相当于100*15的分组交织
end
%*****************************************************************************%
% % % 编码和交织后将二进制序列转换成十进制序列,准备输入到G函数进行映射
% % if (BPH~=1)
% % for i=1:2*N
% % register=interleave_outcome((i-1)*BPH+1:i*BPH);
% % source_coded_ba0(i)=change2deci(register,2); % 其实是把BPH个二进制比特转化为一个2^BPH进制的符号,这里BPH=2,就是把2个二进制比特转化为1个四进制符号(是对的)
% % end
% % else
% % source_coded_ba0=interleave_outcome(1:2*N);% source_coded_ba0 是行矢量,source_coded_ba0是编码输出十进制符号序列,长度为2N个符号(补2L个0符号之前)
% % end
%
% source_coded=[source_coded_ba0,zeros(1,L)];% G函数映射前补上L个0符号
% 编码和交织后将GF(16)域的16进制符号转化为普通域的fanout=4进制的符号,准备输入到G函数进行映射
codedsequence_quart=zeros(1,3000);
for i=1:1500
for j=0:2^M-1
if (codedsequence_GF(i)==j)
codedsequence_quart(2*(i-1)+1:2*i)=quarter_table(j+1,:);
end
end
end % 把1500个编码后的16进制 RS 码元都转化为3000个4进制符号,准备输入到G函数进行映射(通过测试,结果正确,可以转化)
source_coded=[codedsequence_quart,zeros(1,L)]; % G函数映射前补上L=2个0符号,总共3002个4进制符号进行G函数映射
% ******************* G 函数实现部分 ************************ %
% 先定义三个关键矩阵"nextstates" "output" "input"
nextstate=zeros(number_of_states,fanout); % nextstate矩阵:行代表网格图中的各状态(一一对应),列与输入移位寄存器的信息符号一一对应,
% 矩阵中存储的内容是与当前状态和输入符号对应的下一状态号(即存储网格图的状态转移规则)
output=zeros(number_of_states,fanout); % output矩阵:行代表网格图中的各状态(一一对应),列与输入移位寄存器的信息符号一一对应,
% 矩阵中存储的内容是与当前状态和输入符号对应的网格图分支转移输出(分支转移输出是跳频频率号)
input=zeros(number_of_states,number_of_states);
number_of_out=number_of_states*fanout;% 跳频频点数Nt
for i=0:number_of_states-1
for j=0:fanout-1
[next_state,out_put]=G_func1(i,j,L,fanout);
nextstate(i+1,j+1)=next_state;
output(i+1,j+1)=out_put;
input(i+1,next_state+1)=j;
end
end
% ********************************************************* %
depth_of_trellis=length(source_coded);%*******************
Eb=1;
% Es=Eb*BPH*(1/2);
Es=Eb*(14/15);% Es=Eb*BPH*(k/n)
diversi=2; % 2重分集
Ec=Es/diversi;% deversi阶分集后每个码片的能量
% Ej0=(Eb*number_of_out*Q)/(Eb_to_Nj);% 每个多音干扰的能量Ej0
Ej0=(Es*number_of_out)/(BPH*Q*Eb_to_Nj); % 每跳时间内每个多音干扰的能量Ej0
sgma=sqrt(Eb/(2*Eb_to_No));% AWGN的均方根
demod_input=zeros(number_of_out,depth_of_trellis);
f=zeros(1,depth_of_trellis);
% rc=zeros(1,number_of_out); % 移到下面去了
% rs=zeros(1,number_of_out);
D=0; % D 记录网格图的当前状态,这里初始状态是0状态
% % *******************信道和非相干解调部分:加多音干扰和噪声,然后非相干解调****************** %
% for i=1:depth_of_trellis % i表示网格图的时间走势
% f(i)=output(D+1,source_coded(i)+1); % f(i)是i时刻的分支转移输出,即i时刻的跳频频率号,频率号范围是[0,number_of_out-1]而不是[1,number_of_out]
% thyta=2*pi*rand;% 干扰音与跳频信号的相对相位
% % J=randint(1,Q,number_of_out); % J 矩阵中存放Q个干扰音所在的频率号,干扰音所在频率号范围也是[0,number_of_out-1]而不是[1,number_of_out]
% J=gen_multijammer(Q,number_of_out);% J 矩阵中存放Q个干扰音所在的频率号,干扰音所在频率号范围也是[0,number_of_out-1]而不是[1,number_of_out]
% for j=0:number_of_out-1
% if (j==f(i))
% rc(j+1)=sqrt(Es)+sgma*randn;
% rs(j+1)=sgma*randn;
% else
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -