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