IT++ Logo Newcom Logo

turbo.cpp

Go to the documentation of this file.
00001 
00033 #include <itpp/comm/turbo.h>
00034 
00035 
00036 namespace itpp { 
00037 
00038   void Turbo_Codec::set_parameters(ivec gen1, ivec gen2, int constraint_length, const ivec &interleaver_sequence,
00039                                    int in_iterations, std::string in_metric, double in_logmax_scale_factor, 
00040                                    bool in_adaptive_stop,  LLR_calc_unit in_llrcalc)
00041   {
00042     //Set the input parameters:
00043     iterations          = in_iterations;
00044     interleaver_size    = interleaver_sequence.size();
00045     Nuncoded            = interleaver_size;
00046     logmax_scale_factor = in_logmax_scale_factor;
00047     adaptive_stop       = in_adaptive_stop;
00048 
00049     //Check the decoding metric
00050     if (in_metric=="LOGMAX") {
00051       metric = "LOGMAX";
00052     } else if (in_metric=="LOGMAP") {
00053       metric = "LOGMAP";
00054     } else if (in_metric=="MAP") {
00055       metric = "MAP";
00056     } else if (in_metric=="TABLE") {
00057       metric = "TABLE";
00058     } else {
00059       it_error("Turbo_Codec::set_parameters: The decoder metric must be either MAP, LOGMAP or LOGMAX");
00060     }
00061 
00062     if (logmax_scale_factor != 1.0) {
00063       it_assert(metric=="LOGMAX","Turbo_Codec::set_parameters: logmax_scale_factor can only be used together with LOGMAX decoding");
00064     }
00065 
00066     //The RSC Encoders:
00067     rscc1.set_generator_polynomials(gen1, constraint_length);
00068     rscc2.set_generator_polynomials(gen2, constraint_length);
00069     n1 = gen1.length()-1; //Number of parity bits from rscc1
00070     n2 = gen2.length()-1; //Number of parity bits from rscc2
00071     n_tot = 1 + n1 + n2;  //Total number of parity bits and systematic bits
00072 
00073     //Set the number of tail bits:
00074     m_tail = constraint_length - 1;
00075 
00076     //Calculate the number of coded bits per code-block:
00077     Ncoded = Nuncoded * n_tot + m_tail * (1 + n1) + m_tail * (1 + n2); 
00078 
00079     //Set the interleaver sequence
00080     bit_interleaver.set_interleaver_depth(interleaver_size);
00081     float_interleaver.set_interleaver_depth(interleaver_size);
00082     bit_interleaver.set_interleaver_sequence(interleaver_sequence);
00083     float_interleaver.set_interleaver_sequence(interleaver_sequence);
00084 
00085     //Default value of the channel reliability scaling factor is 1
00086     Lc = 1.0;
00087 
00088     // LLR algebra table   
00089     rscc1.set_llrcalc(in_llrcalc);
00090     rscc2.set_llrcalc(in_llrcalc);
00091 
00092   }
00093 
00094   void Turbo_Codec::set_interleaver(const ivec &interleaver_sequence)
00095   {
00096     interleaver_size = interleaver_sequence.size();
00097     Nuncoded = interleaver_size;
00098 
00099     //Calculate the number of coded bits per code-block:
00100     Ncoded = Nuncoded * n_tot + m_tail * (1 + n1) + m_tail * (1 + n2); 
00101   
00102     //Set the interleaver sequence
00103     bit_interleaver.set_interleaver_depth(interleaver_size);
00104     float_interleaver.set_interleaver_depth(interleaver_size);
00105     bit_interleaver.set_interleaver_sequence(interleaver_sequence);
00106     float_interleaver.set_interleaver_sequence(interleaver_sequence);
00107   }
00108 
00109   void Turbo_Codec::set_metric(std::string in_metric, double in_logmax_scale_factor, LLR_calc_unit in_llrcalc)
00110   {
00111     logmax_scale_factor = in_logmax_scale_factor;
00112   
00113     //Check the decoding metric
00114     if (in_metric=="LOGMAX") {
00115       metric = "LOGMAX";
00116     } else if (in_metric=="LOGMAP") {
00117       metric = "LOGMAP";
00118     } else if (in_metric=="MAP") {
00119       metric = "MAP";
00120     } else if (in_metric=="TABLE") {
00121       metric = "TABLE";
00122     } else {
00123       it_error("Turbo_Codec::set_metric: The decoder metric must be either MAP, LOGMAP or LOGMAX");
00124     }
00125 
00126     rscc1.set_llrcalc(in_llrcalc);
00127     rscc2.set_llrcalc(in_llrcalc);
00128   }
00129 
00130   void Turbo_Codec::set_iterations(int in_iterations)
00131   {
00132     iterations = in_iterations;
00133   }
00134 
00135   void Turbo_Codec::set_adaptive_stop(bool in_adaptive_stop)
00136   {
00137     adaptive_stop = in_adaptive_stop;
00138   }
00139 
00140   void Turbo_Codec::set_awgn_channel_parameters(double in_Ec, double in_N0)
00141   {
00142     Ec = in_Ec;
00143     N0 = in_N0;
00144     Lc = 4.0*sqrt(Ec)/N0;
00145   }
00146 
00147   void Turbo_Codec::set_scaling_factor(double in_Lc)
00148   {
00149     Lc = in_Lc;
00150   }
00151 
00152 
00153   void Turbo_Codec::encode(const bvec &input, bvec &output) 
00154   {
00155     //Local variables:
00156     int i, k, j, no_blocks;
00157     long count;
00158     bvec input_bits, in1, in2, tail1, tail2, out;
00159     bmat parity1, parity2;
00160   
00161     //Initializations:
00162     no_blocks = input.length() / Nuncoded;
00163     output.set_size(no_blocks*Ncoded,false);
00164   
00165     //Set the bit counter to zero:
00166     count = 0;
00167 
00168     //Encode all code blocks:
00169     for (i=0; i<no_blocks; i++) {
00170     
00171       //Encode one block
00172       input_bits = input.mid(i*Nuncoded,Nuncoded);
00173       encode_block(input_bits, in1, in2, parity1, parity2);
00174     
00175       //The data part:
00176       for (k=0; k<Nuncoded; k++) {
00177         output(count) = in1(k); count++;                                //Systematic bits
00178         for (j=0; j<n1; j++) { output(count) = parity1(k,j); count++; } //Parity-1 bits
00179         for (j=0; j<n2; j++) { output(count) = parity2(k,j); count++; } //Parity-2 bits
00180       }
00181     
00182       //The first tail:
00183       for (k=0; k<m_tail; k++) {
00184         output(count) = in1(Nuncoded+k); count++;                                //First systematic tail bit
00185         for (j=0; j<n1; j++) { output(count) = parity1(Nuncoded+k,j); count++; } //Parity-1 tail bits
00186       }
00187 
00188       //The second tail:
00189       for (k=0; k<m_tail; k++) {
00190         output(count) = in2(Nuncoded+k); count++;                                //Second systematic tail bit
00191         for (j=0; j<n2; j++) { output(count) = parity2(Nuncoded+k,j); count++; } //Parity-2 tail bits
00192       }
00193 
00194     } 
00195   
00196   }
00197 
00198   void Turbo_Codec::decode(const vec &received_signal, bvec &decoded_bits, const bvec &true_bits)
00199   {
00200     ivec nrof_used_iterations;
00201     decode(received_signal, decoded_bits, nrof_used_iterations, true_bits);
00202   }
00203 
00204   void Turbo_Codec::decode(const vec &received_signal, bvec &decoded_bits, ivec &nrof_used_iterations, 
00205                            const bvec &true_bits)
00206   {
00207 
00208     if ( (n1==1) && (n2==1) && (metric!="MAP") ) {
00209       //This is a speed optimized decoder for R=1/3 (log domain metrics only)
00210       decode_n3(received_signal, decoded_bits, nrof_used_iterations, true_bits);
00211     } else {
00212 
00213       //Local variables:
00214       vec rec, rec_syst1, rec_syst2;
00215       mat rec_parity1, rec_parity2;
00216       bmat decoded_bits_i;
00217       int no_blocks, i, j, k, nrof_used_iterations_i;
00218       long count;
00219       bool CHECK_TRUE_BITS;
00220     
00221       //Initilaizations:
00222       no_blocks = received_signal.length() / Ncoded;
00223       decoded_bits.set_size(no_blocks * Nuncoded, false);
00224       decoded_bits_i.set_size(iterations, no_blocks * Nuncoded, false);
00225       rec_syst1.set_size(Nuncoded+m_tail,false);
00226       rec_syst2.set_size(Nuncoded+m_tail,false); rec_syst2.clear();
00227       rec_parity1.set_size(Nuncoded+m_tail,n1,false); 
00228       rec_parity2.set_size(Nuncoded+m_tail,n2,false); 
00229       nrof_used_iterations.set_size(no_blocks, false);
00230 
00231       //Check the vector true_bits:
00232       if (true_bits.size()>1) {
00233         it_assert(true_bits.size() == (Nuncoded * no_blocks),"Turbo_Codec::decode: Wrong size of input vectors");
00234         CHECK_TRUE_BITS = true;
00235       } else {
00236         CHECK_TRUE_BITS = false;
00237       }
00238     
00239       //Set the bit counter to zero:
00240       count = 0;
00241     
00242       //Itterate over all received code blocks:
00243       for (i=0; i<no_blocks; i++) {
00244       
00245         //The data part:
00246         for (k=0; k<Nuncoded; k++) {
00247           rec_syst1(k) = received_signal(count); count++;                               //Systematic bit
00248           for (j=0; j<n1; j++) { rec_parity1(k,j) = received_signal(count); count++; }  //Parity-1 bits
00249           for (j=0; j<n2; j++) { rec_parity2(k,j) = received_signal(count); count++; }  //Parity-2 bits
00250         }
00251       
00252         //The first tail:
00253         for (k=0; k<m_tail; k++) {
00254           rec_syst1(Nuncoded+k) = received_signal(count); count++;                               //Tail 1 systematic bit
00255           for (j=0; j<n1; j++) { rec_parity1(Nuncoded+k,j) = received_signal(count); count++; }  //Tail 1 parity-1 bits
00256         }
00257       
00258         //The second tail:
00259         for (k=0; k<m_tail; k++) {
00260           rec_syst2(Nuncoded+k) = received_signal(count); count++;                              //Tail2 systematic bit
00261           for (j=0; j<n2; j++) { rec_parity2(Nuncoded+k,j) = received_signal(count); count++; } //Tali2 parity-2 bits
00262         }
00263       
00264         //Scale the input data if necessary:
00265         if (Lc != 1.0) {
00266           rec_syst1 *= Lc;
00267           rec_syst2 *= Lc;
00268           rec_parity1 *= Lc;
00269           rec_parity2 *= Lc;
00270         }
00271       
00272         //Decode the block:
00273         if (CHECK_TRUE_BITS) {
00274           decode_block(rec_syst1, rec_syst2, rec_parity1, rec_parity2, decoded_bits_i, 
00275                        nrof_used_iterations_i, true_bits.mid(i*Nuncoded,Nuncoded) );
00276           nrof_used_iterations(i) = nrof_used_iterations_i;
00277         } else {
00278           decode_block(rec_syst1, rec_syst2, rec_parity1, rec_parity2, decoded_bits_i, nrof_used_iterations_i);
00279           nrof_used_iterations(i) = nrof_used_iterations_i;
00280         }
00281       
00282         //Put the decoded bits in the output vector:
00283         decoded_bits.replace_mid(i*Nuncoded, decoded_bits_i.get_row(iterations-1)); 
00284       
00285       }
00286 
00287     }
00288 
00289   }
00290 
00291   void Turbo_Codec::encode_block(const bvec &input, bvec &in1, bvec &in2, bmat &parity1, bmat &parity2)
00292   {
00293     //Local variables:
00294     bvec tail1, tail2, interleaved_input;
00295 
00296     //Error check:
00297     it_assert(input.length() == Nuncoded,"Turbo_Codec::encode_block: Parameter error in Nuncoded.");
00298 
00299     //Initializations:
00300     tail1.set_size(m_tail, false);                tail1.clear();
00301     tail2.set_size(m_tail, false);                tail2.clear();
00302     parity1.set_size(Nuncoded+m_tail, n1, false); parity1.clear();
00303     parity2.set_size(Nuncoded+m_tail, n2, false); parity2.clear();
00304     interleaved_input.set_size(Nuncoded, false);  interleaved_input.clear();
00305 
00306     //The first encoder:
00307     rscc1.encode_tail(input, tail1, parity1);
00308 
00309     //The interleaver:
00310     bit_interleaver.interleave(input,interleaved_input); 
00311 
00312     //The second encoder:
00313     rscc2.encode_tail(interleaved_input, tail2, parity2);
00314 
00315     //The input vectors used to the two constituent encoders:
00316     in1 = concat(input,tail1);
00317     in2 = concat(interleaved_input,tail2);
00318   
00319   }
00320 
00321   void Turbo_Codec::decode_block(const vec &rec_syst1, const vec &rec_syst2, const mat &rec_parity1, 
00322                                  const mat &rec_parity2, bmat &decoded_bits_i, int &nrof_used_iterations_i, 
00323                                  const bvec &true_bits)
00324   {
00325     //Local variables:
00326     int i;
00327     long count, l, k;
00328     vec extrinsic_input, extrinsic_output, int_rec_syst1, int_rec_syst, tmp;
00329     vec deint_rec_syst2, rec_syst, sub_rec_syst, Le12, Le21, Le12_int, Le21_int, L, tail1, tail2;
00330     bool CHECK_TRUE_BITS, CONTINUE; 
00331 
00332     //Size initializations:
00333     decoded_bits_i.set_size(iterations, Nuncoded, false); 
00334     Le12.set_size(Nuncoded+m_tail, false);
00335     Le21.set_size(Nuncoded+m_tail, false); Le21.zeros();
00336 
00337     //Calculate the interleaved and the deinterleaved sequences:
00338     float_interleaver.interleave(rec_syst1.left(interleaver_size),int_rec_syst1);
00339     float_interleaver.deinterleave(rec_syst2.left(interleaver_size),deint_rec_syst2);
00340   
00341     //Combine the results from rec_syst1 and rec_syst2 (in case some bits are transmitted several times)
00342     rec_syst = rec_syst1.left(interleaver_size) + deint_rec_syst2;
00343     int_rec_syst = rec_syst2.left(interleaver_size) + int_rec_syst1;
00344   
00345     //Get the two tails
00346     tail1 = rec_syst1.right(m_tail);
00347     tail2 = rec_syst2.right(m_tail);
00348 
00349     //Form the input vectors (including tails) to the two decoders:
00350     rec_syst = concat(rec_syst,tail1);
00351     int_rec_syst = concat(int_rec_syst,tail2);
00352   
00353     // Check the vector true_bits
00354     if (true_bits.size()>1) {
00355       it_assert(true_bits.size()==Nuncoded,"Turbo_Codec::decode_block: Illegal size of input vector true_bits");
00356       CHECK_TRUE_BITS = true;
00357     } else {
00358       CHECK_TRUE_BITS = false;
00359     }
00360 
00361     if (CHECK_TRUE_BITS) {
00362       it_assert(adaptive_stop==false,
00363                 "Turbo_Codec::decode_block: You can not stop iterations both adaptively and on true bits");
00364     }
00365 
00366     // Do the iterative decoding:
00367     nrof_used_iterations_i = iterations;
00368     for (i=0; i<iterations; i++) {
00369 
00370       // Decode Code 1
00371       if (metric=="MAP") {
00372         rscc1.map_decode(rec_syst, rec_parity1, Le21, Le12, true);
00373       } else if ((metric=="LOGMAX") || (metric=="LOGMAP") || (metric=="TABLE")) {
00374         rscc1.log_decode(rec_syst, rec_parity1, Le21, Le12, true, metric);
00375         if (logmax_scale_factor != 1.0) {
00376           Le12 *= logmax_scale_factor;
00377         }
00378       } else {
00379         it_error("Turbo_Codec::decode_block: Illegal metric value");
00380       }
00381 
00382       // Interleave the extrinsic information:
00383       float_interleaver.interleave(Le12.left(interleaver_size),tmp);
00384       Le12_int = concat(tmp,zeros(Le12.size()-interleaver_size));
00385 
00386       // Decode Code 2
00387       if (metric=="MAP") {
00388         rscc2.map_decode(int_rec_syst, rec_parity2, Le12_int, Le21_int, true);
00389       } else if ((metric=="LOGMAX") || (metric=="LOGMAP")  || (metric=="TABLE")) {
00390         rscc2.log_decode(int_rec_syst, rec_parity2, Le12_int, Le21_int, true, metric);
00391         if (logmax_scale_factor != 1.0) {
00392           Le21_int *= logmax_scale_factor;
00393         }
00394       } else {
00395         it_error("Turbo_Codec::decode_block: Illegal metric value");
00396       }
00397 
00398       // De-interleave the extrinsic information:
00399       float_interleaver.deinterleave(Le21_int.left(interleaver_size),tmp);
00400       Le21 = concat(tmp,zeros(Le21_int.size()-interleaver_size));  
00401 
00402       // Take bit decisions
00403       L = rec_syst + Le21 + Le12;
00404       count = 0;
00405       for (l=0; l<Nuncoded; l++) {
00406         (L(l)>0.0) ? (decoded_bits_i(i,count) = bin(0)) : (decoded_bits_i(i,count) = bin(1));
00407         count++;
00408       }
00409  
00410       //Check if it is possible to stop iterating early:
00411       CONTINUE = true;
00412       if (i<(iterations-1)) {
00413       
00414         if (CHECK_TRUE_BITS) {
00415           CONTINUE = false;
00416           for (k=0; k<Nuncoded; k++) { if (true_bits(k) != decoded_bits_i(i,k)) { CONTINUE = true; break; } }
00417         } 
00418       
00419         if ((adaptive_stop) && (i>0)) {
00420           CONTINUE = false;
00421           for (k=0; k<Nuncoded; k++) { if (decoded_bits_i(i-1,k) != decoded_bits_i(i,k)) { CONTINUE = true; break; } }
00422         }
00423         
00424       }
00425       
00426       //Check if iterations shall continue:
00427       if (CONTINUE==false) {
00428         //Copy the results from current iteration to all following iterations:
00429         for (k=(i+1); k<iterations; k++) {
00430           decoded_bits_i.set_row(k, decoded_bits_i.get_row(i) );
00431           nrof_used_iterations_i = i+1;
00432         }
00433         break;
00434       }
00435 
00436     }
00437   
00438   }
00439 
00440   void Turbo_Codec::decode_n3(const vec &received_signal, bvec &decoded_bits, ivec &nrof_used_iterations, 
00441                               const bvec &true_bits)
00442   {
00443     //Local variables:
00444     vec rec, rec_syst1, int_rec_syst1, rec_syst2;
00445     vec rec_parity1, rec_parity2;
00446     vec extrinsic_input, extrinsic_output, Le12, Le21, Le12_int, Le21_int, L;
00447     bvec temp_decoded_bits;
00448     int no_blocks, i, j, k, l, nrof_used_iterations_i;
00449     long count, count_out;
00450     bool CHECK_TRUE_BITS, CONTINUE;
00451 
00452     //Initializations:
00453     no_blocks = received_signal.length() / Ncoded;
00454     decoded_bits.set_size(no_blocks * Nuncoded, false);
00455     rec_syst1.set_size(Nuncoded+m_tail,false);
00456     rec_syst2.set_size(Nuncoded+m_tail,false); rec_syst2.clear();
00457     rec_parity1.set_size(Nuncoded+m_tail,false); 
00458     rec_parity2.set_size(Nuncoded+m_tail,false); 
00459     temp_decoded_bits.set_size(Nuncoded,false);
00460     decoded_bits_previous_iteration.set_size(Nuncoded,false);
00461     nrof_used_iterations.set_size(no_blocks, false);
00462 
00463     //Size initializations:
00464     Le12.set_size(Nuncoded, false);
00465     Le21.set_size(Nuncoded, false);
00466 
00467     //Set the bit counter to zero:
00468     count = 0;
00469     count_out = 0;
00470 
00471     // Check the vector true_bits
00472     if (true_bits.size()>1) {
00473       it_assert(true_bits.size()==Nuncoded*no_blocks,"Turbo_Codec::decode_n3: Illegal size of input vector true_bits");
00474       CHECK_TRUE_BITS = true;
00475     } else {
00476       CHECK_TRUE_BITS = false;
00477     }
00478 
00479     if (CHECK_TRUE_BITS) {
00480       it_assert(adaptive_stop==false,
00481                 "Turbo_Codec::decode_block: You can not stop iterations both adaptively and on true bits");
00482     }
00483 
00484     //Iterate over all received code blocks:
00485     for (i=0; i<no_blocks; i++) {
00486 
00487       //Reset extrinsic data:
00488       Le21.zeros();
00489 
00490       //The data part:
00491       for (k=0; k<Nuncoded; k++) {
00492         rec_syst1(k)   = received_signal(count); count++; //Systematic bit
00493         rec_parity1(k) = received_signal(count); count++; //Parity-1 bits
00494         rec_parity2(k) = received_signal(count); count++; //Parity-2 bits
00495       }
00496 
00497       //The first tail:
00498       for (k=0; k<m_tail; k++) {
00499         rec_syst1(Nuncoded+k)   = received_signal(count); count++; //Tail 1 systematic bit
00500         rec_parity1(Nuncoded+k) = received_signal(count); count++; //Tail 1 parity-1 bits
00501       }
00502 
00503       //The second tail:
00504       for (k=0; k<m_tail; k++) {
00505         rec_syst2(Nuncoded+k)   = received_signal(count); count++; //Tail2 systematic bit
00506         rec_parity2(Nuncoded+k) = received_signal(count); count++; //Tali2 parity-2 bits
00507       }
00508     
00509       float_interleaver.interleave(rec_syst1.left(Nuncoded),int_rec_syst1);
00510       rec_syst2.replace_mid(0,int_rec_syst1);
00511 
00512       //Scale the input data if necessary:
00513       if (Lc != 1.0) {
00514         rec_syst1   *= Lc;
00515         rec_syst2   *= Lc;
00516         rec_parity1 *= Lc;
00517         rec_parity2 *= Lc;
00518       }
00519 
00520       //Decode the block:
00521       CONTINUE = true;
00522       nrof_used_iterations_i = iterations;
00523       for (j=0; j<iterations; j++) {
00524 
00525         rscc1.log_decode_n2(rec_syst1, rec_parity1, Le21, Le12, true, metric);
00526         if (logmax_scale_factor!=1.0) { Le12 *= logmax_scale_factor; }
00527         float_interleaver.interleave(Le12,Le12_int);
00528 
00529         rscc2.log_decode_n2(rec_syst2, rec_parity2, Le12_int, Le21_int, true, metric);
00530         if (logmax_scale_factor!=1.0) { Le21_int *= logmax_scale_factor; }
00531         float_interleaver.deinterleave(Le21_int,Le21);
00532 
00533         if (adaptive_stop) {
00534           L = rec_syst1.left(Nuncoded) + Le21.left(Nuncoded) + Le12.left(Nuncoded);
00535           for (l=0; l<Nuncoded; l++) { (L(l)>0.0) ? (temp_decoded_bits(l) = bin(0)) : (temp_decoded_bits(l) = bin(1)); }
00536           if (j==0) { decoded_bits_previous_iteration = temp_decoded_bits; } else {
00537             if (temp_decoded_bits==decoded_bits_previous_iteration ) {
00538               CONTINUE = false;
00539             } else if (j<(iterations-1)) {
00540               decoded_bits_previous_iteration = temp_decoded_bits;
00541             }
00542           }
00543         }
00544 
00545         if (CHECK_TRUE_BITS) {
00546           L = rec_syst1.left(Nuncoded) + Le21.left(Nuncoded) + Le12.left(Nuncoded);
00547           for (l=0; l<Nuncoded; l++) { (L(l)>0.0) ? (temp_decoded_bits(l) = bin(0)) : (temp_decoded_bits(l) = bin(1)); }
00548           if (temp_decoded_bits == true_bits.mid(i*Nuncoded,Nuncoded)) {
00549             CONTINUE = false;
00550           }
00551         }
00552         
00553         if (CONTINUE==false) { nrof_used_iterations_i = j+1; break; }
00554 
00555       }
00556       
00557       //Take final bit decisions
00558       L = rec_syst1.left(Nuncoded) + Le21.left(Nuncoded) + Le12.left(Nuncoded);
00559       for (l=0; l<Nuncoded; l++) { 
00560         (L(l)>0.0) ? (decoded_bits(count_out) = bin(0)) : (decoded_bits(count_out) = bin(1)); count_out++; 
00561       }
00562     
00563       nrof_used_iterations(i) = nrof_used_iterations_i;
00564 
00565     }
00566 
00567   }
00568 
00569   ivec wcdma_turbo_interleaver_sequence(int interleaver_size)
00570   {
00571     const int MAX_INTERLEAVER_SIZE = 5114;
00572     const int MIN_INTERLEAVER_SIZE = 40;
00573     int K;  //Interleaver size
00574     int R;  //Number of rows of rectangular matrix
00575     int C;  //Number of columns of rectangular matrix
00576     int p;  //Prime number
00577     int v;  //Primitive root
00578     ivec s; //Base sequence for intra-row permutation
00579     ivec q; //Minimum prime integers
00580     ivec r; //Permuted prime integers
00581     ivec T; //Inter-row permutation pattern
00582     imat U; //Intra-row permutation patter
00583     ivec I; //The interleaver sequence
00584     ivec primes, roots, Pat1, Pat2, Pat3, Pat4, Isort;
00585     int i, j, qj, temp, row, col, index, count;
00586   
00587     if (interleaver_size > MAX_INTERLEAVER_SIZE) {
00588 
00589       I = sort_index(randu(interleaver_size));
00590       return I;
00591   
00592     } else {
00593 
00594       p = 0;
00595       v = 0;
00596 
00597       //Check the range of the interleaver size:
00598       it_assert(interleaver_size <= MAX_INTERLEAVER_SIZE,"wcdma_turbo_interleaver_sequence: The interleaver size is to large");
00599       it_assert(interleaver_size >= MIN_INTERLEAVER_SIZE,"wcdma_turbo_interleaver_sequence: The interleaver size is to small");
00600     
00601       K = interleaver_size;
00602     
00603       //Definitions of primes and associated primitive roots:
00604       primes = "2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97 101 103 107 109 113 127 131 137 139 149 151 157 163 167 173 179 181 191 193 197 199 211 223 227 229 233 239 241 251 257";
00605       roots = "0 0 0 3 2 2 3 2 5 2 3 2 6 3 5 2 2 2 2 7 5 3 2 3 5 2 5 2 6 3 3 2 3 2 2 6 5 2 5 2 2 2 19 5 2 3 2 3 2 6 3 7 7 6 3";
00606     
00607       //Determine R
00608       if ((K>=40) && (K<=159)) {
00609         R = 5;
00610       } else if ( ((K>=160)&&(K<=200)) || ((K>=481)&&(K<=530)) ) {
00611         R = 10;
00612       } else {
00613         R = 20;
00614       }
00615     
00616       //Determine C
00617       if ((K>=481)&&(K<=530)) {
00618         p = 53;
00619         v = 2;
00620         C = p;
00621       } else {
00622         //Find minimum prime p such that (p+1) - K/R >= 0 ...
00623         for (i=0; i<primes.length(); i++) {
00624           if ( ( double(primes(i)+1) - double(K)/double(R) ) >= 0.0 ) {
00625             p = primes(i);
00626             v = roots(i);
00627             break;
00628           }
00629         }
00630         //... and etermine C such that
00631         if ( ( double(p) - double(K)/double(R) ) >= 0.0 ) {
00632           if ( ( double(p) - 1.0 - double(K)/double(R) ) >= 0.0 ) {
00633             C = p-1;
00634           } else {
00635             C = p;
00636           }
00637         } else {
00638           C = p+1;
00639         }
00640       }
00641     
00642       //Construct the base sequencs s for intra-row permutaions
00643       s.set_size(p-1,false);
00644       s.clear();
00645       s(0) = 1;
00646       for (i=1; i<=(p-2); i++) {
00647         s(i) = (int) mod((long) (v * s(i-1)), (long) p);
00648       }
00649     
00650       //Let q(0) = 1 be the first prime integer in {q(j)}, and select the consecutive 
00651       //minimum prime integers {q(j)}, j = 1, 2, ..., (R-1) such that gcd( q(j), p-1) == 1, q(j) > 6, and q(j) > q(j-1)
00652       q.set_size(R, false);
00653       q.clear();
00654       q(0) = 1;
00655       for (j=1; j<=(R-1); j++) {
00656         for (i=0; i<primes.length(); i++) {
00657           qj = primes(i);
00658           if ( (qj>6) && (qj>q(j-1)) ) {
00659             if (gcd((long) qj,(long) (p-1)) == 1) {
00660               q(j) = qj;
00661               break;
00662             }
00663           }
00664         }
00665       }
00666     
00667       //Definitions of Pat1, Pat2, Pat3, and Pat4:
00668       Pat1 = "19 9 14 4 0 2 5 7 12 18 10 8 13 17 3 1 16 6 15 11";
00669       Pat2 = "19 9 14 4 0 2 5 7 12 18 16 13 17 15 3 1 6 11 8 10";
00670       Pat3 = "9 8 7 6 5 4 3 2 1 0";
00671       Pat4 = "4 3 2 1 0";
00672     
00673       //T(j) is the inter-row permutation patters defined as one of the following four
00674       //kinds of patterns: Pat1, Pat2, Pat3, and Pat4 depending on the number of input bits K
00675       if (K>=3211) {
00676         T = Pat1;
00677       } else if (K>=3161) {
00678         T = Pat2;
00679       } else if (K>=2481) {
00680         T = Pat1;
00681       } else if (K>=2281) {
00682         T = Pat2;
00683       } else if (K>=531) {
00684         T = Pat1;
00685       } else if (K>=481) {
00686         T = Pat3;
00687       } else if (K>=201) {
00688         T = Pat1;
00689       } else if (K>=160) {
00690         T = Pat3;
00691       } else {
00692         T = Pat4;
00693       }
00694     
00695       //Permute {q(j)} to make {r(j)} such that r(T(j)) = q(j), j = 0, 1, ..., (R-1),
00696       //where T(j) indicates the original row position of the j-th permuted row
00697       r.set_size(R,false);
00698       r.clear();
00699       for (j=0; j<=(R-1); j++) {
00700         r( T(j) ) = q(j);
00701       }
00702     
00703       //U(j,i) is the input bit position of i-th output after the permutation of j-th row
00704       //Perform the j-th (j=0, 1, 2, ..., (R-1)) intra-row permutation as
00705       U.set_size(R,C,false);
00706       U.clear();
00707       if (C==p) {    
00708         for (j=0; j<=(R-1); j++) {
00709           for (i=0; i<=(p-2); i++) {
00710             U(j,i) = s( mod((long) (i*r(j)), (long) (p-1)) );
00711           }
00712           U(j,p-1) = 0;
00713         }
00714       } else if (C==(p+1)) {
00715         for (j=0; j<=(R-1); j++) {
00716           for (i=0; i<=(p-2); i++) {
00717             U(j,i) = s( mod((long) (i*r(j)),(long) (p-1)) ); 
00718           }
00719           U(j,p-1) = 0;
00720           U(j,p) = p;
00721         }
00722         if (K == (C*R)) {
00723           temp = U(R-1,p);
00724           U(R-1,p) = U(R-1,0);
00725           U(R-1,0) = temp;
00726         }
00727       } else if (C==(p-1)) {
00728         for (j=0; j<=(R-1); j++) {
00729           for (i=0; i<=(p-2); i++) {
00730             U(j,i) = s( mod((long) (i*r(j)),(long) (p-1)) ) - 1;
00731           }
00732         }
00733       }
00734     
00735       //Calculate the interleaver sequence:
00736       I.set_size(K, false);
00737       I.clear();
00738       count = 0;
00739       for (i=0; i<C; i++) {
00740         for (j=0; j<R; j++) {
00741           row = T(j);
00742           col = U(row,i);
00743           index = row * C + col;
00744           if (index < K) {
00745             I(count) = index;
00746             count++;
00747           }
00748         }
00749       }
00750     
00751       return I;    
00752     }
00753   }
00754 
00755 } // namespace itpp
SourceForge Logo

Generated on Sat Aug 25 23:40:55 2007 for IT++ by Doxygen 1.5.2