00001 00029 #include <itpp/comm/punct_convcode.h> 00030 00031 00032 namespace itpp 00033 { 00034 00035 // --------------------- Punctured_Conv_Code -------------------------------- 00036 00037 //------- Protected functions ----------------------- 00038 int Punctured_Convolutional_Code::weight(int state, int input, int time) 00039 { 00040 int i, j, temp, out, w = 0, shiftreg = state; 00041 00042 shiftreg = shiftreg | (int(input) << m); 00043 for (j = 0; j < n; j++) { 00044 if (puncture_matrix(j, time) == bin(1)) { 00045 out = 0; 00046 temp = shiftreg & gen_pol(j); 00047 for (i = 0; i < K; i++) { 00048 out ^= (temp & 1); 00049 temp = temp >> 1; 00050 } 00051 w += out; 00052 } 00053 } 00054 return w; 00055 } 00056 00057 int Punctured_Convolutional_Code::weight_reverse(int state, int input, int time) 00058 { 00059 int i, j, temp, out, w = 0, shiftreg = state; 00060 00061 shiftreg = shiftreg | (int(input) << m); 00062 for (j = 0; j < n; j++) { 00063 if (puncture_matrix(j, Period - 1 - time) == bin(1)) { 00064 out = 0; 00065 temp = shiftreg & gen_pol_rev(j); 00066 for (i = 0; i < K; i++) { 00067 out ^= (temp & 1); 00068 temp = temp >> 1; 00069 } 00070 w += out; 00071 } 00072 } 00073 return w; 00074 } 00075 00076 void Punctured_Convolutional_Code::weight(int state, int &w0, int &w1, int time) 00077 { 00078 int i, j, temp, out, shiftreg = state; 00079 w0 = 0; 00080 w1 = 0; 00081 00082 shiftreg = shiftreg | (1 << m); 00083 for (j = 0; j < n; j++) { 00084 if (puncture_matrix(j, time) == bin(1)) { 00085 out = 0; 00086 temp = shiftreg & gen_pol(j); 00087 for (i = 0; i < m; i++) { 00088 out ^= (temp & 1); 00089 temp = temp >> 1; 00090 } 00091 w0 += out; 00092 w1 += out ^(temp & 1); 00093 } 00094 } 00095 } 00096 00097 void Punctured_Convolutional_Code::weight_reverse(int state, int &w0, int &w1, int time) 00098 { 00099 int i, j, temp, out, shiftreg = state; 00100 w0 = 0; 00101 w1 = 0; 00102 00103 shiftreg = shiftreg | (1 << m); 00104 for (j = 0; j < n; j++) { 00105 if (puncture_matrix(j, Period - 1 - time) == bin(1)) { 00106 out = 0; 00107 temp = shiftreg & gen_pol_rev(j); 00108 for (i = 0; i < m; i++) { 00109 out ^= (temp & 1); 00110 temp = temp >> 1; 00111 } 00112 w0 += out; 00113 w1 += out ^(temp & 1); 00114 } 00115 } 00116 } 00117 00118 //------- Public functions ----------------------- 00119 00120 void Punctured_Convolutional_Code::set_puncture_matrix(const bmat &pmatrix) 00121 { 00122 it_error_if((pmatrix.rows() != n) || (pmatrix.cols() == 0), "Wrong size of puncture matrix"); 00123 puncture_matrix = pmatrix; 00124 Period = puncture_matrix.cols(); 00125 00126 int p, j; 00127 total = 0; 00128 00129 for (j = 0; j < n; j++) { 00130 for (p = 0; p < Period; p++) 00131 total = total + (int)(puncture_matrix(j, p)); 00132 } 00133 rate = (double)Period / total; 00134 } 00135 00136 void Punctured_Convolutional_Code::encode(const bvec &input, bvec &output) 00137 { 00138 switch (cc_method) { 00139 case Trunc: 00140 encode_trunc(input, output); 00141 break; 00142 case Tail: 00143 encode_tail(input, output); 00144 break; 00145 case Tailbite: 00146 encode_tailbite(input, output); 00147 break; 00148 default: 00149 encode_tail(input, output); 00150 break; 00151 }; 00152 } 00153 00154 void Punctured_Convolutional_Code::encode_trunc(const bvec &input, bvec &output) 00155 { 00156 Convolutional_Code::encode_trunc(input, output); 00157 00158 int nn = 0, i, p = 0, j; 00159 00160 for (i = 0; i < int(output.size() / n); i++) { 00161 for (j = 0; j < n; j++) { 00162 if (puncture_matrix(j, p) == bin(1)) { 00163 output(nn) = output(i * n + j); 00164 nn++; 00165 } 00166 } 00167 p = (p + 1) % Period; 00168 } 00169 output.set_size(nn, true); 00170 } 00171 00172 void Punctured_Convolutional_Code::encode_tail(const bvec &input, bvec &output) 00173 { 00174 Convolutional_Code::encode_tail(input, output); 00175 00176 int nn = 0, i, p = 0, j; 00177 00178 for (i = 0; i < int(output.size() / n); i++) { 00179 for (j = 0; j < n; j++) { 00180 if (puncture_matrix(j, p) == bin(1)) { 00181 output(nn) = output(i * n + j); 00182 nn++; 00183 } 00184 } 00185 p = (p + 1) % Period; 00186 } 00187 output.set_size(nn, true); 00188 } 00189 00190 void Punctured_Convolutional_Code::encode_tailbite(const bvec &input, bvec &output) 00191 { 00192 Convolutional_Code::encode_tailbite(input, output); 00193 00194 int nn = 0, i, p = 0, j; 00195 00196 for (i = 0; i < int(output.size() / n); i++) { 00197 for (j = 0; j < n; j++) { 00198 if (puncture_matrix(j, p) == bin(1)) { 00199 output(nn) = output(i * n + j); 00200 nn++; 00201 } 00202 } 00203 p = (p + 1) % Period; 00204 } 00205 output.set_size(nn, true); 00206 } 00207 00208 00209 // --------------- Hard-decision decoding is not implemented -------------------------------- 00210 void Punctured_Convolutional_Code::decode(const bvec &, bvec &) 00211 { 00212 it_error("Punctured_Convolutional_Code::decode(bvec, bvec); hard-decision decoding is not implemented"); 00213 } 00214 00215 bvec Punctured_Convolutional_Code::decode(const bvec &) 00216 { 00217 it_error("Punctured_Convolutional_Code::decode(bvec, bvec); hard-decision decoding is not implemented"); 00218 return bvec(); 00219 } 00220 00221 /* 00222 Decode a block of encoded data using specified method 00223 */ 00224 void Punctured_Convolutional_Code::decode(const vec &received_signal, bvec &output) 00225 { 00226 switch (cc_method) { 00227 case Trunc: 00228 decode_trunc(received_signal, output); 00229 break; 00230 case Tail: 00231 decode_tail(received_signal, output); 00232 break; 00233 case Tailbite: 00234 decode_tailbite(received_signal, output); 00235 break; 00236 default: 00237 decode_tail(received_signal, output); 00238 break; 00239 }; 00240 } 00241 00242 00243 // Viterbi decoder using TruncLength (=5*K if not specified) 00244 void Punctured_Convolutional_Code::decode_trunc(const vec &received_signal, bvec &output) 00245 { 00246 int nn = 0, i = 0, p = received_signal.size() / total, j; 00247 00248 int temp_size = p * Period * n; 00249 // number of punctured bits in a fraction of the puncture matrix 00250 // correcponding to the end of the received_signal vector 00251 p = received_signal.size() - p * total; 00252 // precise calculation of the number of unpunctured bits 00253 // (in the above fraction of the puncture matrix) 00254 while (p > 0) { 00255 for (j = 0; j < n; j++) { 00256 if (puncture_matrix(j, nn) == bin(1)) 00257 p--; 00258 } 00259 nn++; 00260 } 00261 temp_size += n * nn; 00262 if (p != 0) { 00263 it_warning("Punctured_Convolutional_Code::decode(): Improper length of " 00264 "the received punctured block, dummy bits have been added"); 00265 } 00266 00267 vec temp(temp_size); 00268 nn = 0; 00269 j = 0; 00270 p = 0; 00271 00272 while (nn < temp.size()) { 00273 if ((puncture_matrix(j, p) == bin(1)) && (i < received_signal.size())) { 00274 temp(nn) = received_signal(i); 00275 i++; 00276 } 00277 else { // insert dummy symbols with the same contribution for 0 and 1 00278 temp(nn) = 0; 00279 } 00280 00281 nn++; 00282 j++; 00283 00284 if (j == n) { 00285 j = 0; 00286 p = (p + 1) % Period; 00287 } 00288 } // while 00289 00290 Convolutional_Code::decode_trunc(temp, output); 00291 } 00292 00293 // Viterbi decoder using TruncLength (=5*K if not specified) 00294 void Punctured_Convolutional_Code::decode_tail(const vec &received_signal, bvec &output) 00295 { 00296 int nn = 0, i = 0, p = received_signal.size() / total, j; 00297 00298 int temp_size = p * Period * n; 00299 // number of punctured bits in a fraction of the puncture matrix 00300 // correcponding to the end of the received_signal vector 00301 p = received_signal.size() - p * total; 00302 // precise calculation of the number of unpunctured bits 00303 // (in the above fraction of the puncture matrix) 00304 while (p > 0) { 00305 for (j = 0; j < n; j++) { 00306 if (puncture_matrix(j, nn) == bin(1)) 00307 p--; 00308 } 00309 nn++; 00310 } 00311 temp_size += n * nn; 00312 if (p != 0) { 00313 it_warning("Punctured_Convolutional_Code::decode_tail(): Improper length " 00314 "of the received punctured block, dummy bits have been added"); 00315 } 00316 00317 vec temp(temp_size); 00318 00319 nn = 0; 00320 j = 0; 00321 p = 0; 00322 00323 while (nn < temp.size()) { 00324 if ((puncture_matrix(j, p) == bin(1)) && (i < received_signal.size())) { 00325 temp(nn) = received_signal(i); 00326 i++; 00327 } 00328 else { // insert dummy symbols with same contribution for 0 and 1 00329 temp(nn) = 0; 00330 } 00331 00332 nn++; 00333 j++; 00334 00335 if (j == n) { 00336 j = 0; 00337 p = (p + 1) % Period; 00338 } 00339 } // while 00340 00341 Convolutional_Code::decode_tail(temp, output); 00342 } 00343 00344 // Decode a block of encoded data where encode_tailbite has been used. Tries all start states. 00345 void Punctured_Convolutional_Code::decode_tailbite(const vec &received_signal, bvec &output) 00346 { 00347 int nn = 0, i = 0, p = received_signal.size() / total, j; 00348 00349 int temp_size = p * Period * n; 00350 // number of punctured bits in a fraction of the puncture matrix 00351 // correcponding to the end of the received_signal vector 00352 p = received_signal.size() - p * total; 00353 // precise calculation of the number of unpunctured bits 00354 // (in the above fraction of the puncture matrix) 00355 while (p > 0) { 00356 for (j = 0; j < n; j++) { 00357 if (puncture_matrix(j, nn) == bin(1)) 00358 p--; 00359 } 00360 nn++; 00361 } 00362 temp_size += n * nn; 00363 if (p != 0) { 00364 it_warning("Punctured_Convolutional_Code::decode_tailbite(): Improper " 00365 "length of the received punctured block, dummy bits have been " 00366 "added"); 00367 } 00368 00369 vec temp(temp_size); 00370 00371 nn = 0; 00372 j = 0; 00373 p = 0; 00374 00375 while (nn < temp.size()) { 00376 if ((puncture_matrix(j, p) == bin(1)) && (i < received_signal.size())) { 00377 temp(nn) = received_signal(i); 00378 i++; 00379 } 00380 else { // insert dummy symbols with same contribution for 0 and 1 00381 temp(nn) = 0; 00382 } 00383 00384 nn++; 00385 j++; 00386 00387 if (j == n) { 00388 j = 0; 00389 p = (p + 1) % Period; 00390 } 00391 } // while 00392 00393 Convolutional_Code::decode_tailbite(temp, output); 00394 } 00395 00396 00397 /* 00398 Calculate the inverse sequence 00399 00400 Assumes that encode_tail is used in the encoding process. Returns false if there is an 00401 error in the coded sequence (not a valid codeword). Does not check that the tail forces 00402 the encoder into the zeroth state. 00403 */ 00404 bool Punctured_Convolutional_Code::inverse_tail(const bvec coded_sequence, bvec &input) 00405 { 00406 int state = 0, zero_state, one_state, zero_temp, one_temp, i, j, p = 0, prev_pos = 0, no_symbols; 00407 int block_length = 0; 00408 bvec zero_output(n), one_output(n), temp_output(n); 00409 00410 block_length = coded_sequence.size() * Period / total - m; 00411 00412 it_error_if(block_length <= 0, "The input sequence is to short"); 00413 input.set_length(block_length, false); // not include the tail in the output 00414 00415 p = 0; 00416 00417 for (i = 0; i < block_length; i++) { 00418 zero_state = state; 00419 one_state = state | (1 << m); 00420 no_symbols = 0; 00421 for (j = 0; j < n; j++) { 00422 if (puncture_matrix(j, p) == bin(1)) { 00423 zero_temp = zero_state & gen_pol(j); 00424 one_temp = one_state & gen_pol(j); 00425 zero_output(no_symbols) = xor_int_table(zero_temp); 00426 one_output(no_symbols) = xor_int_table(one_temp); 00427 no_symbols++; 00428 } 00429 } 00430 if (coded_sequence.mid(prev_pos, no_symbols) == zero_output.left(no_symbols)) { 00431 input(i) = bin(0); 00432 state = zero_state >> 1; 00433 } 00434 else if (coded_sequence.mid(prev_pos, no_symbols) == one_output.left(no_symbols)) { 00435 input(i) = bin(1); 00436 state = one_state >> 1; 00437 } 00438 else 00439 return false; 00440 00441 prev_pos += no_symbols; 00442 p = (p + 1) % Period; 00443 } 00444 00445 return true; 00446 } 00447 00448 00449 00450 00451 /* 00452 It is possible to do this more efficiently by registrating all (inputs,states,Periods) 00453 that has zero metric (just and with the generators). After that build paths from a input=1 state. 00454 */ 00455 bool Punctured_Convolutional_Code::catastrophic(void) 00456 { 00457 int max_stack_size = 50000; 00458 ivec S_stack(max_stack_size), t_stack(max_stack_size); 00459 int start, S, W0, W1, S0, S1, pos, t = 0; 00460 int stack_pos = -1; 00461 00462 for (pos = 0; pos < Period; pos++) { // test all start positions 00463 for (start = 0; start < (1 << m); start++) { 00464 stack_pos = -1; 00465 S = start; 00466 t = 0; 00467 00468 node1: 00469 if (t > (1 << m) * Period) { return true; } 00470 S0 = next_state(S, 0); 00471 S1 = next_state(S, 1); 00472 weight(S, W0, W1, (pos + t) % Period); 00473 if (W1 > 0) { goto node0; } 00474 if ((W0 == 0) && (S0 == start) && (((pos + t + 1) % Period) == pos)) { return true; } 00475 if ((W0 == 0) && (S0 == 0) && (S != 0)) { return true; } 00476 if ((S0 != 0) && (W0 == 0)) { 00477 stack_pos++; 00478 if (stack_pos >= max_stack_size) { 00479 max_stack_size = 2 * max_stack_size; 00480 S_stack.set_size(max_stack_size, true); 00481 t_stack.set_size(max_stack_size, true); 00482 } 00483 S_stack(stack_pos) = S0; 00484 t_stack(stack_pos) = t + 1; 00485 } 00486 if ((W1 == 0) && (S1 == start) && (((pos + t + 1) % Period) == pos)) { return true; } 00487 S = S1; 00488 t++; 00489 goto node1; 00490 00491 node0: 00492 if (W0 > 0) goto stack; 00493 if ((W0 == 0) && (S0 == start) && (((pos + t + 1) % Period) == pos)) { return true; } 00494 if ((W0 == 0) && (S0 == 0) && (S != 0)) { return true; } 00495 if (S0 != 0) { 00496 S = S0; 00497 t++; 00498 goto node1; 00499 } 00500 else { 00501 goto stack; 00502 } 00503 00504 stack: 00505 if (stack_pos >= 0) { 00506 S = S_stack(stack_pos); 00507 t = t_stack(stack_pos); 00508 stack_pos--; 00509 goto node1; 00510 } 00511 } 00512 } 00513 return false; 00514 } 00515 00516 void Punctured_Convolutional_Code::distance_profile(ivec &dist_prof, int start_time, int dmax, bool reverse) 00517 { 00518 int max_stack_size = 50000; 00519 ivec S_stack(max_stack_size), W_stack(max_stack_size), t_stack(max_stack_size); 00520 00521 int stack_pos = -1, t, S, W, W0, w0, w1; 00522 00523 00524 dist_prof.zeros(); 00525 dist_prof += dmax; // just select a big number! 00526 if (reverse) 00527 W = weight_reverse(0, 1, start_time); 00528 else 00529 W = weight(0, 1, start_time); 00530 00531 S = next_state(0, 1); // start from zero state and a one input 00532 t = 0; 00533 dist_prof(0) = W; 00534 00535 node1: 00536 if (reverse) 00537 weight_reverse(S, w0, w1, (start_time + t + 1) % Period); 00538 else 00539 weight(S, w0, w1, (start_time + t + 1) % Period); 00540 00541 if (t < m) { 00542 W0 = W + w0; 00543 if (W0 < dist_prof(m)) { // store node0 (S, W0, and t+1) 00544 stack_pos++; 00545 if (stack_pos >= max_stack_size) { 00546 max_stack_size = 2 * max_stack_size; 00547 S_stack.set_size(max_stack_size, true); 00548 W_stack.set_size(max_stack_size, true); 00549 t_stack.set_size(max_stack_size, true); 00550 } 00551 S_stack(stack_pos) = next_state(S, 0); 00552 W_stack(stack_pos) = W0; 00553 t_stack(stack_pos) = t + 1; 00554 } 00555 } 00556 else goto stack; 00557 00558 W += w1; 00559 if (W > dist_prof(m)) 00560 goto stack; 00561 00562 t++; 00563 S = next_state(S, 1); 00564 00565 if (W < dist_prof(t)) 00566 dist_prof(t) = W; 00567 00568 if (t == m) goto stack; 00569 goto node1; 00570 00571 00572 stack: 00573 if (stack_pos >= 0) { 00574 // get root node from stack 00575 S = S_stack(stack_pos); 00576 W = W_stack(stack_pos); 00577 t = t_stack(stack_pos); 00578 stack_pos--; 00579 00580 if (W < dist_prof(t)) 00581 dist_prof(t) = W; 00582 00583 if (t == m) goto stack; 00584 00585 goto node1; 00586 } 00587 } 00588 00589 int Punctured_Convolutional_Code::fast(Array<ivec> &spectrum, int start_time, int dfree, int no_terms, 00590 int d_best_so_far, bool test_catastrophic) 00591 { 00592 int cat_treshold = (1 << m) * Period; 00593 int i, j, t = 0; 00594 ivec dist_prof(K), dist_prof_rev(K), dist_prof_temp(K); 00595 00596 //calculate distance profile 00597 distance_profile(dist_prof, start_time, dfree); 00598 distance_profile(dist_prof_rev, 0, dfree, true); // for the reverse code 00599 00600 int dist_prof_rev0 = dist_prof_rev(0); 00601 00602 bool reverse = true; // true = use reverse dist_prof 00603 00604 // choose the lowest dist_prof over all periods 00605 for (i = 1; i < Period; i++) { 00606 distance_profile(dist_prof_temp, i, dfree, true); 00607 // switch if new is lower 00608 for (j = 0; j < K; j++) { 00609 if (dist_prof_temp(j) < dist_prof_rev(j)) { 00610 dist_prof_rev(j) = dist_prof_temp(j); 00611 } 00612 } 00613 } 00614 00615 dist_prof_rev0 = dist_prof(0); 00616 dist_prof = dist_prof_rev; 00617 00618 int d = dfree + no_terms - 1; 00619 int max_stack_size = 50000; 00620 ivec S_stack(max_stack_size), W_stack(max_stack_size), b_stack(max_stack_size); 00621 ivec c_stack(max_stack_size), t_stack(max_stack_size); 00622 int stack_pos = -1; 00623 00624 // F1. 00625 int mf = 1, b = 1; 00626 spectrum.set_size(2); 00627 spectrum(0).set_size(dfree + no_terms, 0); // ad 00628 spectrum(1).set_size(dfree + no_terms, 0); // cd 00629 spectrum(0).zeros(); 00630 spectrum(1).zeros(); 00631 int S, S0, S1, W0, W1, w0, w1, c = 0; 00632 S = next_state(0, 1); // start in zero state and one input 00633 int W = d - dist_prof_rev0; 00634 t = 1; 00635 00636 F2: 00637 S0 = next_state(S, 0); 00638 S1 = next_state(S, 1); 00639 00640 if (reverse) { 00641 weight(S, w0, w1, (start_time + t) % Period); 00642 } 00643 else { 00644 weight_reverse(S, w0, w1, (start_time + t) % Period); 00645 } 00646 W0 = W - w0; 00647 W1 = W - w1; 00648 00649 if (mf < m) goto F6; 00650 00651 //F3: 00652 if (W0 >= 0) { 00653 spectrum(0)(d - W0)++; 00654 spectrum(1)(d - W0) += b; 00655 } 00656 //Test on d_best_so_far 00657 if ((d - W0) < d_best_so_far) { return -1; } 00658 00659 F4: 00660 if ((W1 < dist_prof(m - 1)) || (W < dist_prof(m))) goto F5; 00661 // select node 1 00662 if (test_catastrophic && (W == W1)) { 00663 c++; 00664 if (c > cat_treshold) 00665 return 0; 00666 } 00667 else { 00668 c = 0; 00669 } 00670 00671 W = W1; 00672 S = S1; 00673 t++; 00674 mf = 1; 00675 b++; 00676 goto F2; 00677 00678 F5: 00679 if (stack_pos == -1) goto end; 00680 // get node from stack 00681 S = S_stack(stack_pos); 00682 W = W_stack(stack_pos); 00683 b = b_stack(stack_pos); 00684 c = c_stack(stack_pos); 00685 t = t_stack(stack_pos); 00686 stack_pos--; 00687 mf = 1; 00688 goto F2; 00689 00690 F6: 00691 if (W0 < dist_prof(m - mf - 1)) goto F4; 00692 00693 //F7: 00694 if ((W1 >= dist_prof(m - 1)) && (W >= dist_prof(m))) { 00695 // save node 1 00696 if (test_catastrophic && (stack_pos > 40000)) 00697 return 0; 00698 00699 stack_pos++; 00700 if (stack_pos >= max_stack_size) { 00701 max_stack_size = 2 * max_stack_size; 00702 S_stack.set_size(max_stack_size, true); 00703 W_stack.set_size(max_stack_size, true); 00704 b_stack.set_size(max_stack_size, true); 00705 c_stack.set_size(max_stack_size, true); 00706 t_stack.set_size(max_stack_size, true); 00707 } 00708 S_stack(stack_pos) = S1; 00709 W_stack(stack_pos) = W1; 00710 b_stack(stack_pos) = b + 1; // because gone to node 1 00711 c_stack(stack_pos) = c; 00712 t_stack(stack_pos) = t + 1; 00713 } 00714 // select node 0 00715 S = S0; 00716 t++; 00717 if (test_catastrophic && (W == W0)) { 00718 c++; 00719 if (c > cat_treshold) 00720 return false; 00721 } 00722 else { 00723 c = 0; 00724 } 00725 00726 W = W0; 00727 mf++; 00728 goto F2; 00729 00730 end: 00731 return 1; 00732 } 00733 00734 void Punctured_Convolutional_Code::calculate_spectrum(Array<ivec> &spectrum, int dmax, int no_terms) 00735 { 00736 Array<ivec> temp_spectra(2); 00737 spectrum.set_size(2); 00738 spectrum(0).set_size(dmax + no_terms, false); 00739 spectrum(1).set_size(dmax + no_terms, false); 00740 spectrum(0).zeros(); 00741 spectrum(1).zeros(); 00742 00743 for (int pos = 0; pos < Period; pos++) { 00744 calculate_spectrum(temp_spectra, pos, dmax, no_terms); 00745 spectrum(0) += temp_spectra(0); 00746 spectrum(1) += temp_spectra(1); 00747 } 00748 } 00749 00750 void Punctured_Convolutional_Code::calculate_spectrum(Array<ivec> &spectrum, int time, int dmax, int no_terms, int block_length) 00751 { 00752 imat Ad_states(1 << (K - 1), dmax + no_terms), Cd_states(1 << m, dmax + no_terms); 00753 imat Ad_temp(1 << (K - 1), dmax + no_terms), Cd_temp(1 << m, dmax + no_terms); 00754 ivec mindist(1 << (K - 1)), mindist_temp(1 << m); 00755 00756 spectrum.set_size(2); 00757 spectrum(0).set_size(dmax + no_terms, false); 00758 spectrum(1).set_size(dmax + no_terms, false); 00759 spectrum(0).zeros(); 00760 spectrum(1).zeros(); 00761 Ad_states.zeros(); 00762 Cd_states.zeros(); 00763 mindist.zeros(); 00764 int wmax = dmax + no_terms; 00765 ivec visited_states(1 << m), visited_states_temp(1 << m); 00766 bool proceede, expand_s1; 00767 int t, d, w0, w1, s, s0, s1 = 0, s_start; 00768 00769 // initial phase. Evolv trellis up to time K. 00770 visited_states.zeros(); // 0 = false 00771 00772 s1 = next_state(0, 1); // Start in state 0 and 1 input 00773 visited_states(s1) = 1; // 1 = true. 00774 w1 = weight(0, 1, time); 00775 Ad_states(s1, w1) = 1; 00776 Cd_states(s1, w1) = 1; 00777 mindist(s1) = w1; 00778 00779 if (block_length > 0) { 00780 s0 = next_state(0, 0); 00781 visited_states(s0) = 1; // 1 = true. Expand also the zero-path 00782 w0 = weight(0, 0, time); 00783 Ad_states(s0, w0) = 1; 00784 Cd_states(s0, w0) = 0; 00785 mindist(s0) = w0; 00786 s_start = 0; 00787 } 00788 else { 00789 s_start = 1; 00790 } 00791 00792 t = 1; 00793 proceede = true; 00794 while (proceede) { 00795 proceede = false; 00796 Ad_temp.zeros(); 00797 Cd_temp.zeros(); 00798 mindist_temp.zeros(); 00799 visited_states_temp.zeros(); //false 00800 00801 for (s = s_start; s < (1 << m); s++) { 00802 00803 if (visited_states(s) == 1) { 00804 if ((mindist(s) >= 0) && (mindist(s) < wmax)) { 00805 proceede = true; 00806 weight(s, w0, w1, (time + t) % Period); 00807 00808 s0 = next_state(s, 0); 00809 for (d = mindist(s); d < (wmax - w0); d++) { 00810 Ad_temp(s0, d + w0) += Ad_states(s, d); 00811 Cd_temp(s0, d + w0) += Cd_states(s, d); 00812 } 00813 if (visited_states_temp(s0) == 0) { mindist_temp(s0) = mindist(s) + w0; } 00814 else { mindist_temp(s0) = ((mindist(s) + w0) < mindist_temp(s0)) ? mindist(s) + w0 : mindist_temp(s0); } 00815 visited_states_temp(s0) = 1; //true 00816 00817 expand_s1 = false; 00818 if ((block_length == 0) || (t < (block_length - (K - 1)))) { 00819 expand_s1 = true; 00820 } 00821 00822 if (expand_s1) { 00823 s1 = next_state(s, 1); 00824 for (d = mindist(s); d < (wmax - w1); d++) { 00825 Ad_temp(s1, d + w1) += Ad_states(s, d); 00826 Cd_temp(s1, d + w1) += Cd_states(s, d) + Ad_states(s, d); 00827 } 00828 if (visited_states_temp(s1) == 0) { mindist_temp(s1) = mindist(s) + w1; } 00829 else { mindist_temp(s1) = ((mindist(s) + w1) < mindist_temp(s1)) ? mindist(s) + w1 : mindist_temp(s1); } 00830 visited_states_temp(s1) = 1; //true 00831 } 00832 } 00833 } 00834 } 00835 00836 Ad_states = Ad_temp; 00837 Cd_states = Cd_temp; 00838 if (block_length == 0) { 00839 spectrum(0) += Ad_temp.get_row(0); 00840 spectrum(1) += Cd_temp.get_row(0); 00841 } 00842 visited_states = visited_states_temp; 00843 mindist = elem_mult(mindist_temp, visited_states); 00844 t++; 00845 if ((t > block_length) && (block_length > 0)) { proceede = false; } 00846 } 00847 00848 if (block_length > 0) { 00849 spectrum(0) = Ad_states.get_row(0); 00850 spectrum(1) = Cd_states.get_row(0); 00851 } 00852 00853 } 00854 00855 } // namespace itpp
Generated on Wed Jul 27 2011 16:27:05 for IT++ by Doxygen 1.7.4