Viterbi译码程序代码

1、译码主要部分#include stdafx.h/#define DEBUGvoid deci2bin(in t d, int size, int * b);int bin2deci(in t *b, int size);int nxt_stat(int current_state, int input, int *memory_contents);void init_quantizer(void);void init_adaptive_quant(float es_ovr_nO);int soft_quant(float channel_symbol);int soft_metr ic(in

2、t data, int guess);int quantizer_table256;void sdvd(intg2K, float es_ovr_nO, long channel_length, float*channel_output_vector, int *decoder_output_matrix)int i, j, l, ll; /循环控制变量long t;/时间int memory_contentsK; /记录输入内容int inputTW0T0THEMTW0T0THEM; /对当前状态以及下一个状态映射int outputTW0T0THEM2; /卷积码编码输出矩阵int nex

3、tstateTW0T0THEM2; /下一个状态矩阵int accum_err_metricTW0T0THEM2; /误差累计矩阵int state_historyTW0T0THEMK * 5 + 1; /历史状态表int state_sequenceK * 5 + 1; /状态序列int *channel_output_matrix; /信道输出序列int binary_output 2;int branch_output2; /0 或者 1 的输出分支int m, n, number_of_s tat es, dep th_o f_t rellis, st ep, branch_me tr

4、 ic,sh_p tr, sh_col, x, xx, h, hh, nex t_state, las t_stop;n = 2;/1/2为卷积码传输数据的码率m = K - 1;/寄存器个数number_of_states = (int)pow(2.0, m);/状态个数 number of states = 2(K - 1) = 2mdepth_of_trellis = K * 5;for (i = 0; i number_of_states; i+)for (j = 0; j number_of_states; j+)inpu t ij = 0;/输入数组初始化for (j = 0; j

5、 n; j+)nex tstat eij = 0;/下一个状态数组初始化out pu tij = 0;/输出数组初始化for (j = 0; j = depth_of_trellis; j+)state_historyij = 0;/历史状态数组初始化 state_history416accum_err_metrici0 = 0;/误差累计矩阵第一列初始化为0accum_err_me trici1 = MAXINT;/误差累计矩阵第二列初始化为一个很大的数/*前向纠错简称 FEC(Forward Error Correction),其原理是：发送方将要发送的数据附加上一定的冗余纠错码一并发送，

6、接收方则根据纠错码对数据进行差错检测，如发现差错，由接收方进行纠正*/*产生状态转移矩阵、输出矩阵、输入矩阵*/输入矩阵表示的是FEC编码传输给下一个状态/下一个状态由输入和当前状态给出/输出矩阵for (j = 0; j number_of_states; j+)for (l = 0; l n; l+)nex t_state = nxt_stat(j, l, memory_con tent s);inp utj nex t_state = l;/*计算给定的卷积编码器输出当前状态数和输入值*/branch_output0 = 0;branch_output1 = 0;for (i = 0;

7、i K; i+)branch_ output0 = branch_ output0 memory_con tent si & g0i; branch_ output1 = branch_ output1 memory_con tent si & g1i;nex tstatejl = nex t_state ;/下一个状态out pu tjl = bin2deci(branch_ou tput, 2);/输出十进制#ifdef DEBUGprintf(nlnput:);for (j = 0; j number_of_states; j+)prin tf (n);for (l = 0; l num

8、ber_of_states; l+)printf(%2d , inputjl);printf(n0utput:);for (j = 0; j number_of_states; j+)prin tf (n);for (l = 0; l n; l+)printf(%2d , outputjl);printf (nNext State:);for (j = 0; j number_of_states; j+)prin tf (n);for (l = 0; l n; l+)prin tf (%2d , nex tstatejl);#endifchannel.output_matrix =(int *

9、)malloc(channel_length * sizeof(int);if (channel.output_matrix = NULL)printf(allocation is failure!n); exit(1);prin tf (n);/*信道输出为n行，2列，每行对应于一个通道符号给定的位和每一列对应于一个已编码的位*/ channel_length = channel_length / n;ini t_adap ti ve_quan t( es_ovr_n0);/进行优化，匹配过信噪比的量化矩阵/量化信道输出，将浮点型数据转化为整形for (t = 0; t (channel_l

10、eng th * n); t = t + n)for (i = 0; i n; i+)*(channel.output_matrix + (t / n) + (i * channel_length)= soft_quant(*(channel.output_vector + (t + i);/printf(%d ,*(channel_output_matrix + (t / n) + (i * channel_length);/*结束设置：利用网格遍历开始译码通道，在编码完成后结束*/for (t = 0; t channel_length - m; t+)if (t二 m)/假设从零，所以只

11、是计算路径的所有零状态step = (int)pow(2.0, m - t * 1);/如果不写成2.0，会出现函数重载不明确的错误 elsestep = 1;/利用state_history矩阵作为循环缓冲区sh_ptr = (int)(t + 1) % (depth_of_trellis + 1);/sh_ptr 为 state history 矩阵的指for (j = 0; j number_of_states; j += step) /重复每个可能的卷积编码器的输出组for (l = 0; l 1; binary_output1 = outputjl & 0 x00000001;bra

12、nch_metrie = branch_metrie +abs(*(channel.output_matrix + (0 * channel_length + t) - 7 * binary_output0)+ abs(*( channel.output_matrix + (1 * channel_length + t) - 7 * binary_output1);/选择累加误差最小的if (accum_e rr_metr icnex tstatejl1accum_e rr_metr icj0 + branch_metric)accum_e rr_metr icnex tstatejl1 =

13、accum_e rr_metr icj0 +branch_metric;state _his tor ynex tstatej lsh_p tr = j;/prin tf (s tate _his tor y%d%d=%dn, nex tstatejl, sh_p tr, state _his tor ynex tstatej lsh_p tr); /循环l结束 /j结束，更新网格/accum_err_metrie矩阵第二列移到第一列，第二列标志为一个很大的数 for (j = 0; j number_of_states; j+)accum_e rr_metr icj0 = accum_e r

14、r_metr icj1;/prin tf (accum_e rr_metr ic%d0=%dn, j, accum_e rr_metr icj0); accum_e rr_metr icj1 = MAXINT;/如果网格填充完成，现在需要追踪for (j = 0; j = depth_of_trellis; j+) /初始化状态序列矩阵 state_sequencej = 0;/for找到的最小累积state_history元素/forMAXINT;(j = 0; j (number_of_s tates / 2); j+)if(accum_err_metricj0 accum_err_met

15、ricnumber_of_states if(accum_err_metricj0 accum_err_metricnumber_of_states 1-j0)xx =accum_e rr_metr icj0;xx =accum_e rr_metr icj0;hh =j;hh =j;elsexx =hh =accum_err_metricnumber_of_states elsexx =hh =accum_err_metricnumber_of_states - 1 - j0; number_of_states 1 j;if(xx x)if(xx 0; j一)sh_col = j + (sh_

16、ptr - depth_of_trellis);if (sh_col 0)sh_col = sh_col + depth_of_trellis + 1;state_sequencej - 1 = state_historystate_sequencejsh_col;/找出输入序列对应的状态序列在最佳路径*(decoder_output_matrix + t - depth_of_trellis + 1)=inp utstate _sequenceOs tate _sequence1;/printf(译码输出:%dn, *(decoder_output_matrix + t depth_of_t

17、rellis + 1); /if状态 /结束t循环/译码信道中的数据for (t = channel_length m; t channel_length; t+)sh_p tr = (in t)( t + 1) % (dep th_of_ trellis + 1);las t_s top = number_of_s tates / pow(2.0, t - channel_leng th + m);/不需要考虑输入的 状态是1,所以确定最高可能的状态数是0for (j = 0; j last_stop; j+)branch_metrie = 0;deci2bin(o utputj0, n,

18、binary _output);for (ll = 0; ll n; ll+)branch_metrie = branch_metrie + soft_metric(*(channel_output_matrix + (ll * channel_length + t), binary_outputll);if (accum_e rr_metr icnex tstatej01accum_e rr_metr icj0 + branch_metric)accum_e rr_metr icnex tstatej01 = accum_e rr_metr icj0 + branch_metric;stat

19、e _his tor ynex tstatej 0sh_p tr = j;for (j = 0; j number_of_states; j+)accum_e rr_metr icjO = accum_e rr_metr icj1; accum_e rr_metr icj1 = MAXINT;/对所选路径进行选择for (j = 0; j = depth_of_trellis; j+) state_sequencej = 0;x = accum_e rr_metr ic00;h = 0;for (j = 1; j last_stop; j+)if (accum_e rr_metr icj0 0

20、; j一)sh_col = j + (sh_ptr - depth_of_trellis);if (sh_col 0)sh_col = sh_col + depth_of_trellis + 1;state_sequencej - 1 = state_historystate_sequencejsh_col; *( decoder_output_matrix + t - depth_of_trellis + 1)=inp utstate _sequence0s tate _sequence1; /if条件状态 /结束t循环for (i = 1; i depth_of_trellis - m;

21、i+)*(decoder_output_matrix + channel_length - depth_of_trellis + i)=inputstate_sequenceistate_sequencei + 1; free(channel_output_matrix);return;/初始化三位软判决量化编码器/加入噪声后的量化void init _adap ti ve_quan tf loa t es_ovr_nO)int i, d;floa t es, sn_ra tio, sigma;es = 1;sn_ra tio = (floa t)pow(10.0, (es_ovr_nO /

22、10.0);sigma = (floa t)sq rt(es / (2.0 * sn_ra tio);d = (in t)(32 * 0.5 * sigma);for (i = -128; i (-3 * d); i+)quantizer_tablei + 128 = 7;for (i = (-3 * d); i (-2 * d); i+)quantizer_tablei + 128 = 6;for (i = (-2 * d); i (-1 * d); i+)quantizer_tablei + 128 = 5;for (i = (-1 * d); i 0; i+)quantizer_tabl

23、ei + 128 = 4;for (i = 0; i (1 * d); i+)quantizer_tablei + 128 = 3;for (i = (1 * d); i (2 * d); i+)quantizer_tablei + 128 = 2;for (i = (2 * d); i (3 * d); i+)quantizer_tablei + 128 = 1;for (i = (3 * d); i 128; i+)quantizer_tablei + 128 = 0;/channel_symbol信道中的值为-1或+1的加噪信号int soft_ quan tf loa t channe

24、l_symbol)int x;x = (int)(32.0 * channel_symbol); /则 x 的平均值为-32 或+32if (x 127) x = 127;/大于 127 则输出 127ret urn(quan tizer_ tablex + 128);/查找量化表/* this metric is based on the algorithm given in Michelson and Levesque, page 323. */int soft_metr ic(in t data, int guess)retu rn(abs(da ta - (guess * 7);/当给

25、出当前状态以及输入时，计算下一个状态，并且计算卷积码编码内容int nxt_stat(int current_state, int input, int *memory_contents)int binary_stateK 1; /二进制当前状态int next_state_binaryK 1; /二进制下一个状态int next_state; /十进制当前状态int i;deci2bin(curren t_state, K 1, binary_s tat e);next_state_binaryO = input;for (i = 1; i K 1; i+)next_state_binary

26、i = binary_statei - 1;next_state = bin2deci(next_state_binary, K 1);memory_contents0 = input;for (i = 1; i K; i+)memory_contentsi = binary_statei 1; /编码内容 ret urn (nex t_s tat e);/返回十进制的下一个状态/将十进制数据转化为特殊的二进制，例如：十进制“ 100011 ”输出二进制数据100011 void deci2bin(in t d, int size, int *b) int i;for (i = 0; i =

27、0; i一) d = d 1;bi = d & 0 x01;/将2进制数据转化为特殊的十进制，例如：二进制“ 100011 ”输出十进制数据100011 int bin2deci(in t *b, int size) int i, d;d = 0;for (i = 0; i size; i+)d += bi (size i 1);return (d);程序辅助模块/ viterbi.cpp :定义控制台应用程序的入口点。/#include stdafx.h int* gen_da ta(long da ta_len)/使用指针函数返回生成的输入数组int *array;unsigned lon

28、g t;/系统时间t,产生随机数序列数字srand(unsigned) time(NULL);/头文件是#include stdlib.harray = (int *)malloc(data_len * sizeof (int);/开辟动态数组，为了防止栈溢出 /将产生的随机数序列写入一个输出数组for (t = 0; t 0.5);/*prin tf (%d , *( array + t);if (t+1)%50=0)prin tf (n);*/return array;int* conv(in t g2K, long inp ut_ len, int *i n_array) long t,

29、 tt;long t, tt;int j, k;int *unencoded_data;int shift_regK;int sr_head;int p, q;/指针型未编码数组/移位寄存器/指向一位寄存器的首位/异或门输出变量int *out_array;int m = K - 1;/编码存储长度，也就是寄存器的个数（K为编码约束长度）unencoded_data = (int *)malloc(input_len + m) * sizeof (int);/给将要编码的序列扩展长 度并分配动态内存out _array = (in t *)malloc(2 * (inp ut _len + m

30、) *sizeof (in t);/一个输入对应两个输出，所 以其长度为2 * (input_len + m)if (unencoded_data = NULL) printf(allocation is faliure ); exit(1);for (t = 0; t input_len; t+) /将系统生成的长度为 Input_len 的数据写入 unencoded_data *( unencoded_da ta + t) = *( in_array + t);/printf (%dn, *(unencoded_data + t);for (t = 0; t m; t+)/在 unenc

31、oded_data 需要编码存储个 0*(unencoded_data + input_len + t) = 0;/printf (%dn, *(unencoded_data + input_len+ t);for (j = 0; j K; j+)/初始化移位寄存器shift_regj = 0;sr_head = 0; /移位寄存器的首指针 tt = 0;/初始化通道标志位输出指针/卷积码编码for (t = 0; t input_len + m; t+) shift_regsr_head = * (unencoded_data + t);/移位寄存器首指针指向 unencoded_data

32、p = 0;q = 0;for (j = 0; j K; j+) k = (j + sr_head) % K;p =p shift_regk & gOj;q = qshi ft_ regk & g1j;/将异或门的两个输出端写入out_array *(out _array + tt) = p;tt = tt + 1;*(out _array + tt) = q;tt = tt + 1;sr_head = sr_headT;/使移位寄存器循环if (sr_head 0)sr_head = m;/*for (int i = 0; i = 1);/第三步：若s大于1,则返回第一步if/否则，计算x1

33、(!s) x1 = 0;elsex1 = v1*sqrt(-2.0*log(s) / s); /log(s)默认对数的底为 ereturn (mean + sigma*x1); float* addnoise(float es_ovr_nO, long channel_len, int *in_array)float *out_array;long t;floa t mean, es, sn_ra tio, sigma, signal;mean = 0;es = 1;/固定比特能量为es=1sn_ratio = (float)pow(10, (es_ovr_nO / 10);/信噪比，由 dB

34、 形式转换成普通形式 sigma = (floa t)sqr t(es / (2 * sn_ra tio); /标准差out_array = (float *)malloc(channel_len *sizeof (float);/给通过 BPSK 信道输出 out_array 开辟空间for (t = 0; t channel_len; t+)signal = 1-2 * (* (in_array + t);/如果输入的是0,经过BPSK变为+1 ；如果输入的是1, 经过BPSK变为-1*(out _array + t) = signal + gngauss(mean, sigma); /加

35、入高斯白噪声/prin tf (%fn, *(out _array + t);/*for (int i = 0; i channel_len; i+)prin tf (%f , *(out _array + i);if (i + 1) % 10= 0)prin tf (n);*/return out_array; extern void sdvd(int g2K, float es_ovr_nO, long channel_length, float *channel.output_vector, int *decoder_output_matrix);int _tm ain(in t arg

36、c, _TCHAR * argv)/*long data_len=LENGTH;int a = 2;/给in_array, out_array两个动态数组提供起始地址float b = 1.0;int *in_array ;int *out_array;float *bpsk_array;bpsk_array = &b;in_array=out_array = &a;int g2K = 1, 1, 1 , 1, 0, 1 ;/当 K=3 时，两个子生成元float mean, sigma;printf(数据的长度：dn, data_len);in_array=gen_data(data_len

37、);/* prin tf (%d n, dat a_len);for (int t = 0; t data_len; t+)prin tf (%d , *( in_array + t);*/*out _array=conv(g, dat a_len,in_array);/* for (in t i = 0; i2 *(da ta_len+KT); i+)prin tf (%d , *(out _array + i);*/* mean = 0.0;sigma = 0.1;gngauss(mean, sigma);prin tf (%fn, gngauss(mean, sigma) );*/*lo

38、ng message = LENGTH;long channel_len = 2*(m essage+K-1);/ 信道长度long es_ovr_n0=L0ESN0;bpsk_array=addnoise(es_ovr_n0, channel_len, out_array);for (int i = 0; i channel_len; i+)prin tf (%f , *( bpsk_array + i);if (i % 2 != 0)prin tf (n);*/long iter, t, msg_length, channel_length;int *original;int *encod

39、ed;/* original, encoded, & decoded data arrays */int *decoded;double start,f inish;float * noisy;/噪声序列int m;/ m = K - 1float es_ovr_nO;int number_errors_encoded;#if K = 3int g2K = 1, 1, 1 ,/* 7 */ 1, 0, 1 ;/* 5 */#endif#if K = 5int g2K = 1, 1, 1, 0, 1 ,/* 35 */ 1, 0, 0, 1, 1 ; /* 23 */#endif#if K =

40、3printf(g1 = %d%d%d, g00, g01, g02);printf(ng2 = %d%d%dn, g10, g11, g12);#endif#if K = 5printf(g1 = %d%d %d%d%d, g00, g01, g02, g03, g04);printf(ng2 = n%d%d %d%d%dn, g10, g11, g12, g13, g14);#endifm = K - 1;msg_length = MSG_LEN;channel_length = (msg_length + m) * 2;prin tf(信息 leng th=%dn, msg_leng th);original = (in