Matlab实现 维特比译码
作者:互联网
利用matlab实现维特比译码
话不多说,直接展示代码
- 首先是利用卷积码进行编码
- 然后加入awgn噪声
- 最后利用维特比译码
- 绘制误比特率图像
卷积码进行编码
%% 初始化参数
N = 1000; %序列长度
code_in = randi(2,1,N)-1;
%code_in = [1,1,0,1,1,1,0,0,1];
N = length(code_in);
g1 = [1,1,1];
g2 = [1,0,1];
%% 计算 卷积码
%x_g1
x0 = g1(1) * code_in;
x1 = g1(2) * code_in;
x2 = g1(3)* code_in;
x = xor( x0(2:N), x1(1:N-1) );
x = [ x0(1), x, x1(N) ];
x3 = xor(x(3:N+1),x2(1:N-1));
x_g1 = [x(1:2) ,x3, x2(N) ];
%x_g2
x0 = g2(1) * code_in;
x1 = g2(2) * code_in;
x2 = g2(3)* code_in;
x = xor( x0(2:N), x1(1:N-1) );
x = [ x0(1), x, x1(N) ];
x3 = xor(x(3:N+1),x2(1:N-1));
x_g2 = [x(1:2) ,x3, x2(N) ];
%% 合并两路
x = zeros(1,size(x_g1, 2)+size(x_g2, 2));
x(1:2:end) = x_g1;
x(2:2:end) = x_g2;
x = x(1:length(x) - 4);
维特比译码
function m = viterbi_hard(x)
%%
% a 00
% b 10
% c 01
% d 11
%x = [ 1 1 0 1 0 1 0 0 0 1];
sizex = size(x);
s = sizex(2)/2;
x = [0,0,0,0, x];
% to record the value
val_a = 0;
val_b = 0;
val_c = 0;
val_d = 0;
% graph aa0 ab1 bc0 bd1 ca0 cb1 dc0 dd1
gra = [ 0,0; 1,1; 1,0; 0,1; 1,1; 0,0; 0,1; 1,0 ];
% to record route
ma = zeros(1,s+2);
mb = zeros(1,s+2);
mc = zeros(1,s+2);
md = zeros(1,s+2);
%% stage 1 and 2
val_a = val_a + dis(gra(1,:), x(1:2));
ma(1)=0;
val_b = val_b + dis(gra(2,:), x(1:2));
mb(1)=1;
mc = mb;
md = mb;
val_c = val_b + dis(gra(3,:), x(3:4));
mc(2)=0;
val_d = val_b + dis(gra(4,:), x(3:4));
md(2)=1;
val_a = val_a + dis(gra(1,:), x(3:4));
ma(2)=0;
val_b = val_b + dis(gra(2,:), x(3:4));
mb(2)=1;
for i = 3:s+2
% val_a_t =val_a;
% val_b_t =val_b;
% val_c_t =val_c;
% val_d_t =val_d;
% tempa = ma;
% tempb = mb;
% tempc = mc;
% tempd = md;
% for val_a
if val_a + dis(gra(1,:), x(2*i-1:2*i)) >= val_c + dis(gra(5,:),x(2*i-1:2*i))
tempa = mc;
val_a_t = val_c + dis(gra(5,:),x(2*i-1:2*i));
tempa(i)=0;
else
val_a_t = val_a + dis(gra(1,:),x(2*i-1:2*i));
tempa = ma;
tempa(i)=0;
end
%for val_b
if val_a + dis(gra(2,:), x(2*i-1:2*i)) >= val_c + dis(gra(6,:),x(2*i-1:2*i))
tempb = mc;
val_b_t = val_c + dis(gra(6,:),x(2*i-1:2*i));
tempb(i)=1;
else
val_b_t = val_a + dis(gra(2,:),x(2*i-1:2*i));
tempb = ma;
tempb(i)=1;
end
%for val_c
if val_b + dis(gra(3,:), x(2*i-1:2*i)) >= val_d + dis(gra(7,:),x(2*i-1:2*i))
tempc = md;
val_c_t = val_d + dis(gra(7,:),x(2*i-1:2*i));
tempc(i)=0;
else
val_c_t = val_b + dis(gra(3,:),x(2*i-1:2*i));
tempc = mb;
tempc(i)=0;
end
%for val_c
if val_b + dis(gra(4,:), x(2*i-1:2*i)) >= val_d + dis(gra(8,:),x(2*i-1:2*i))
tempd = md;
val_d_t = val_d + dis(gra(8,:),x(2*i-1:2*i));
tempd(i)=1;
else
val_d_t = val_b + dis(gra(4,:),x(2*i-1:2*i));
tempd = mb;
tempd(i)=1;
end
val_a =val_a_t;
val_b =val_b_t;
val_c =val_c_t;
val_d =val_d_t;
ma = tempa;
mb = tempb;
mc = tempc;
md = tempd;
end
if val_a <= val_b
m = ma;
t = val_a;
else
m = mb;
t = val_a;
end
if val_c <= t
m = mc;
t =val_c;
end
if val_d <= t
m = md;
t = val_d;
end
%% 去掉最开始的00
m = m(3:s+2);
译码中间用到dis 函数,测量两段序列的海明距离
%% function to compute distances between x,y in the form of [0,0] [1,1]
function d = dis(x, y)
d = sum(xor(x,y));
加入噪声,绘制误比特率随snr变化的曲线
errbit = zeros(1,21);
for j = -5:15
y = awgn(x, j, 'measured');
for i = 1:2*N
if y(i) >=0.5
y(i)=1;
else
y(i)=0;
end
end
% 译码
m=viterbi_hard(y);
errbit(j+6) = sum(m~=code_in)/N ;
end
logerr = 10*log10(errbit);
plot( -5:15, logerr );
标签:code,end,val,维特,mb,译码,Matlab,gra,dis 来源: https://blog.csdn.net/fms2000/article/details/117848705