1 /* 2 * Copyright (c) 2007 - 2019 Joseph Gaeddert 3 * 4 * Permission is hereby granted, free of charge, to any person obtaining a copy 5 * of this software and associated documentation files (the "Software"), to deal 6 * in the Software without restriction, including without limitation the rights 7 * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell 8 * copies of the Software, and to permit persons to whom the Software is 9 * furnished to do so, subject to the following conditions: 10 * 11 * The above copyright notice and this permission notice shall be included in 12 * all copies or substantial portions of the Software. 13 * 14 * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR 15 * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, 16 * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE 17 * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER 18 * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, 19 * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN 20 * THE SOFTWARE. 21 */ 22 23 // 24 // liquid.internal.h 25 // 26 // Internal header file for liquid DSP for SDR 27 // 28 // This file includes function declarations which are intended 29 // for internal use 30 // 31 32 #ifndef __LIQUID_INTERNAL_H__ 33 #define __LIQUID_INTERNAL_H__ 34 35 // Configuration file 36 #include "config.h" 37 38 #include <complex.h> 39 #include "liquid.h" 40 41 #if defined HAVE_FEC_H && defined HAVE_LIBFEC 42 # define LIBFEC_ENABLED 1 43 #endif 44 45 46 // 47 // Debugging macros 48 // 49 #define DEBUG_PRINTF_FLOAT(F,STR,I,V) \ 50 fprintf(F,"%s(%4u) = %12.4e;\n",STR,I+1,V) 51 #define DEBUG_PRINTF_CFLOAT(F,STR,I,V) \ 52 fprintf(F,"%s(%4u) = %12.4e +j*%12.4e;\n",STR,I+1,crealf(V),cimagf(V)) 53 54 #define PRINTVAL_FLOAT(X,F) printf(#F,crealf(X)); 55 #define PRINTVAL_CFLOAT(X,F) printf(#F "+j*" #F, crealf(X), cimagf(X)); 56 57 // 58 // MODULE : agc 59 // 60 61 62 // 63 // MODULE : audio 64 // 65 66 67 // 68 // MODULE : buffer 69 // 70 71 72 // 73 // MODULE : dotprod 74 // 75 76 77 // 78 // MODULE : fec (forward error-correction) 79 // 80 81 // checksum / cyclic redundancy check (crc) 82 83 #define CRC8_POLY 0x07 84 #define CRC16_POLY 0x8005 85 #define CRC24_POLY 0x5D6DCB 86 #define CRC32_POLY 0x04C11DB7 87 88 unsigned int checksum_generate_key(unsigned char * _msg, unsigned int _msg_len); 89 unsigned int crc8_generate_key(unsigned char * _msg, unsigned int _msg_len); 90 unsigned int crc16_generate_key(unsigned char * _msg, unsigned int _msg_len); 91 unsigned int crc24_generate_key(unsigned char * _msg, unsigned int _msg_len); 92 unsigned int crc32_generate_key(unsigned char * _msg, unsigned int _msg_len); 93 94 95 // fec : basic object 96 struct fec_s { 97 // common 98 fec_scheme scheme; 99 //unsigned int dec_msg_len; 100 //unsigned int enc_msg_len; 101 float rate; 102 103 // lengths: convolutional, Reed-Solomon 104 unsigned int num_dec_bytes; 105 unsigned int num_enc_bytes; 106 107 // convolutional : internal memory structure 108 unsigned char * enc_bits; 109 void * vp; // decoder object 110 int * poly; // polynomial 111 unsigned int R; // primitive rate, inverted (e.g. R=3 for 1/3) 112 unsigned int K; // constraint length 113 unsigned int P; // puncturing rate (e.g. p=3 for 3/4) 114 int * puncturing_matrix; 115 116 // viterbi decoder function pointers 117 void*(*create_viterbi)(int); 118 //void (*set_viterbi_polynomial)(int*); 119 int (*init_viterbi)(void*,int); 120 int (*update_viterbi_blk)(void*,unsigned char*,int); 121 int (*chainback_viterbi)(void*,unsigned char*,unsigned int,unsigned int); 122 void (*delete_viterbi)(void*); 123 124 // Reed-Solomon 125 int symsize; // symbol size (bits per symbol) 126 int genpoly; // generator polynomial 127 int fcs; // 128 int prim; // 129 int nroots; // number of roots in the polynomial 130 //int ntrials; // 131 unsigned int rspad; // number of implicit padded symbols 132 int nn; // 2^symsize - 1 133 int kk; // nn - nroots 134 void * rs; // Reed-Solomon internal object 135 136 // Reed-Solomon decoder 137 unsigned int num_blocks; // number of blocks: ceil(dec_msg_len / nn) 138 unsigned int dec_block_len; // number of decoded bytes per block: 139 unsigned int enc_block_len; // number of encoded bytes per block: 140 unsigned int res_block_len; // residual bytes in last block 141 unsigned int pad; // padding for each block 142 unsigned char * tblock; // decoder input sequence [size: 1 x n] 143 int * errlocs; // error locations [size: 1 x n] 144 int * derrlocs; // decoded error locations [size: 1 x n] 145 int erasures; // number of erasures 146 147 // encode function pointer 148 void (*encode_func)(fec _q, 149 unsigned int _dec_msg_len, 150 unsigned char * _msg_dec, 151 unsigned char * _msg_enc); 152 153 // decode function pointer 154 void (*decode_func)(fec _q, 155 unsigned int _dec_msg_len, 156 unsigned char * _msg_enc, 157 unsigned char * _msg_dec); 158 159 // decode function pointer (soft decision) 160 void (*decode_soft_func)(fec _q, 161 unsigned int _dec_msg_len, 162 unsigned char * _msg_enc, 163 unsigned char * _msg_dec); 164 }; 165 166 // simple type testing 167 int fec_scheme_is_convolutional(fec_scheme _scheme); 168 int fec_scheme_is_punctured(fec_scheme _scheme); 169 int fec_scheme_is_reedsolomon(fec_scheme _scheme); 170 int fec_scheme_is_hamming(fec_scheme _scheme); 171 int fec_scheme_is_repeat(fec_scheme _scheme); 172 173 // Pass 174 fec fec_pass_create(void *_opts); 175 void fec_pass_destroy(fec _q); 176 void fec_pass_print(fec _q); 177 void fec_pass_encode(fec _q, 178 unsigned int _dec_msg_len, 179 unsigned char * _msg_dec, 180 unsigned char * _msg_enc); 181 void fec_pass_decode(fec _q, 182 unsigned int _dec_msg_len, 183 unsigned char * _msg_enc, 184 unsigned char * _msg_dec); 185 186 // Repeat (3) 187 fec fec_rep3_create(void *_opts); 188 void fec_rep3_destroy(fec _q); 189 void fec_rep3_print(fec _q); 190 void fec_rep3_encode(fec _q, 191 unsigned int _dec_msg_len, 192 unsigned char * _msg_dec, 193 unsigned char * _msg_enc); 194 void fec_rep3_decode(fec _q, 195 unsigned int _dec_msg_len, 196 unsigned char * _msg_enc, 197 unsigned char * _msg_dec); 198 void fec_rep3_decode_soft(fec _q, 199 unsigned int _dec_msg_len, 200 unsigned char * _msg_enc, 201 unsigned char * _msg_dec); 202 203 // Repeat (5) 204 fec fec_rep5_create(void *_opts); 205 void fec_rep5_destroy(fec _q); 206 void fec_rep5_print(fec _q); 207 void fec_rep5_encode(fec _q, 208 unsigned int _dec_msg_len, 209 unsigned char * _msg_dec, 210 unsigned char * _msg_enc); 211 void fec_rep5_decode(fec _q, 212 unsigned int _dec_msg_len, 213 unsigned char * _msg_enc, 214 unsigned char * _msg_dec); 215 void fec_rep5_decode_soft(fec _q, 216 unsigned int _dec_msg_len, 217 unsigned char * _msg_enc, 218 unsigned char * _msg_dec); 219 220 // Hamming(7,4) 221 extern unsigned char hamming74_enc_gentab[16]; 222 extern unsigned char hamming74_dec_gentab[128]; 223 fec fec_hamming74_create(void *_opts); 224 void fec_hamming74_destroy(fec _q); 225 void fec_hamming74_print(fec _q); 226 void fec_hamming74_encode(fec _q, 227 unsigned int _dec_msg_len, 228 unsigned char * _msg_dec, 229 unsigned char * _msg_enc); 230 void fec_hamming74_decode(fec _q, 231 unsigned int _dec_msg_len, 232 unsigned char * _msg_enc, 233 unsigned char * _msg_dec); 234 void fec_hamming74_decode_soft(fec _q, 235 unsigned int _dec_msg_len, 236 unsigned char * _msg_enc, 237 unsigned char * _msg_dec); 238 // soft decoding of one symbol 239 unsigned char fecsoft_hamming74_decode(unsigned char * _soft_bits); 240 241 // Hamming(8,4) 242 extern unsigned char hamming84_enc_gentab[16]; 243 extern unsigned char hamming84_dec_gentab[256]; 244 fec fec_hamming84_create(void *_opts); 245 void fec_hamming84_destroy(fec _q); 246 void fec_hamming84_print(fec _q); 247 void fec_hamming84_encode(fec _q, 248 unsigned int _dec_msg_len, 249 unsigned char * _msg_dec, 250 unsigned char * _msg_enc); 251 void fec_hamming84_decode(fec _q, 252 unsigned int _dec_msg_len, 253 unsigned char * _msg_enc, 254 unsigned char * _msg_dec); 255 void fec_hamming84_decode_soft(fec _q, 256 unsigned int _dec_msg_len, 257 unsigned char * _msg_enc, 258 unsigned char * _msg_dec); 259 // soft decoding of one symbol 260 unsigned char fecsoft_hamming84_decode(unsigned char * _soft_bits); 261 262 // Hamming(12,8) 263 264 unsigned int fec_hamming128_encode_symbol(unsigned int _sym_dec); 265 unsigned int fec_hamming128_decode_symbol(unsigned int _sym_enc); 266 extern unsigned short int hamming128_enc_gentab[256]; // encoding table 267 268 fec fec_hamming128_create(void *_opts); 269 void fec_hamming128_destroy(fec _q); 270 void fec_hamming128_print(fec _q); 271 void fec_hamming128_encode(fec _q, 272 unsigned int _dec_msg_len, 273 unsigned char * _msg_dec, 274 unsigned char * _msg_enc); 275 void fec_hamming128_decode(fec _q, 276 unsigned int _dec_msg_len, 277 unsigned char * _msg_enc, 278 unsigned char * _msg_dec); 279 void fec_hamming128_decode_soft(fec _q, 280 unsigned int _dec_msg_len, 281 unsigned char * _msg_enc, 282 unsigned char * _msg_dec); 283 // soft decoding of one symbol 284 unsigned int fecsoft_hamming128_decode(unsigned char * _soft_bits); 285 extern unsigned char fecsoft_hamming128_n3[256][17]; 286 unsigned int fecsoft_hamming128_decode_n3(unsigned char * _soft_bits); 287 288 289 // Hamming(15,11) 290 unsigned int fec_hamming1511_encode_symbol(unsigned int _sym_dec); 291 unsigned int fec_hamming1511_decode_symbol(unsigned int _sym_enc); 292 293 // Hamming(31,26) 294 unsigned int fec_hamming3126_encode_symbol(unsigned int _sym_dec); 295 unsigned int fec_hamming3126_decode_symbol(unsigned int _sym_enc); 296 297 298 // Golay(24,12) 299 300 unsigned int fec_golay2412_encode_symbol(unsigned int _sym_dec); 301 unsigned int fec_golay2412_decode_symbol(unsigned int _sym_enc); 302 extern unsigned int golay2412_P[12]; 303 extern unsigned int golay2412_Gt[24]; 304 extern unsigned int golay2412_H[12]; 305 306 // multiply input vector with matrix 307 unsigned int golay2412_matrix_mul(unsigned int _v, 308 unsigned int * _A, 309 unsigned int _n); 310 311 // search for p[i] such that w(v+p[i]) <= 2, return -1 on fail 312 int golay2412_parity_search(unsigned int _v); 313 314 fec fec_golay2412_create(void *_opts); 315 void fec_golay2412_destroy(fec _q); 316 void fec_golay2412_print(fec _q); 317 void fec_golay2412_encode(fec _q, 318 unsigned int _dec_msg_len, 319 unsigned char * _msg_dec, 320 unsigned char * _msg_enc); 321 void fec_golay2412_decode(fec _q, 322 unsigned int _dec_msg_len, 323 unsigned char * _msg_enc, 324 unsigned char * _msg_dec); 325 326 // SEC-DED (22,16) 327 328 // compute parity on 16-bit input 329 unsigned char fec_secded2216_compute_parity(unsigned char * _m); 330 331 // compute syndrome on 22-bit input 332 unsigned char fec_secded2216_compute_syndrome(unsigned char * _v); 333 334 // encode symbol 335 // _sym_dec : decoded symbol [size: 2 x 1] 336 // _sym_enc : encoded symbol [size: 3 x 1], _sym_enc[0] has only 6 bits 337 void fec_secded2216_encode_symbol(unsigned char * _sym_dec, 338 unsigned char * _sym_enc); 339 340 // decode symbol, returning 0/1/2 for zero/one/multiple errors detected 341 // _sym_enc : encoded symbol [size: 3 x 1], _sym_enc[0] has only 6 bits 342 // _sym_dec : decoded symbol [size: 2 x 1] 343 int fec_secded2216_decode_symbol(unsigned char * _sym_enc, 344 unsigned char * _sym_dec); 345 346 // estimate error vector, returning 0/1/2 for zero/one/multiple errors detected 347 // _sym_enc : encoded symbol [size: 3 x 1], _sym_enc[0] has only 6 bits 348 // _e_hat : estimated error vector [size: 3 x 1] 349 int fec_secded2216_estimate_ehat(unsigned char * _sym_enc, 350 unsigned char * _e_hat); 351 352 // parity matrix [6 x 16 bits], [6 x 2 bytes] 353 extern unsigned char secded2216_P[12]; 354 355 // syndrome vectors of errors with weight exactly equal to 1 356 extern unsigned char secded2216_syndrome_w1[22]; 357 358 fec fec_secded2216_create(void *_opts); 359 void fec_secded2216_destroy(fec _q); 360 void fec_secded2216_print(fec _q); 361 void fec_secded2216_encode(fec _q, 362 unsigned int _dec_msg_len, 363 unsigned char * _msg_dec, 364 unsigned char * _msg_enc); 365 void fec_secded2216_decode(fec _q, 366 unsigned int _dec_msg_len, 367 unsigned char * _msg_enc, 368 unsigned char * _msg_dec); 369 370 // SEC-DED (39,32) 371 372 // compute parity on 32-bit input 373 unsigned char fec_secded3932_compute_parity(unsigned char * _m); 374 375 // compute syndrome on 39-bit input 376 unsigned char fec_secded3932_compute_syndrome(unsigned char * _v); 377 378 // encode symbol 379 // _sym_dec : decoded symbol [size: 4 x 1] 380 // _sym_enc : encoded symbol [size: 5 x 1], _sym_enc[0] has only 7 bits 381 void fec_secded3932_encode_symbol(unsigned char * _sym_dec, 382 unsigned char * _sym_enc); 383 384 // estimate error vector, returning 0/1/2 for zero/one/multiple errors detected 385 // _sym_enc : encoded symbol [size: 5 x 1], _sym_enc[0] has only 7 bits 386 // _e_hat : estimated error vector [size: 5 x 1] 387 int fec_secded3932_estimate_ehat(unsigned char * _sym_enc, 388 unsigned char * _e_hat); 389 390 // decode symbol, returning 0/1/2 for zero/one/multiple errors detected 391 // _sym_enc : encoded symbol [size: 5 x 1], _sym_enc[0] has only 7 bits 392 // _sym_dec : decoded symbol [size: 4 x 1] 393 int fec_secded3932_decode_symbol(unsigned char * _sym_enc, 394 unsigned char * _sym_dec); 395 396 // parity matrix [7 x 32 bits], [7 x 4 bytes] 397 extern unsigned char secded3932_P[28]; 398 399 // syndrome vectors of errors with weight exactly equal to 1 400 extern unsigned char secded3932_syndrome_w1[39]; 401 402 fec fec_secded3932_create(void *_opts); 403 void fec_secded3932_destroy(fec _q); 404 void fec_secded3932_print(fec _q); 405 void fec_secded3932_encode(fec _q, 406 unsigned int _dec_msg_len, 407 unsigned char * _msg_dec, 408 unsigned char * _msg_enc); 409 void fec_secded3932_decode(fec _q, 410 unsigned int _dec_msg_len, 411 unsigned char * _msg_enc, 412 unsigned char * _msg_dec); 413 414 // SEC-DED (72,64) 415 416 // compute parity byte on 64-byte input 417 unsigned char fec_secded7264_compute_parity(unsigned char * _v); 418 419 // compute syndrome on 72-bit input 420 unsigned char fec_secded7264_compute_syndrome(unsigned char * _v); 421 422 // encode symbol 423 // _sym_dec : input symbol [size: 8 x 1] 424 // _sym_enc : input symbol [size: 9 x 1] 425 void fec_secded7264_encode_symbol(unsigned char * _sym_dec, 426 unsigned char * _sym_enc); 427 428 // estimate error vector, returning 0/1/2 for zero/one/multiple errors detected 429 // _sym_enc : encoded symbol [size: 9 x 1] 430 // _e_hat : estimated error vector [size: 9 x 1] 431 int fec_secded7264_estimate_ehat(unsigned char * _sym_enc, 432 unsigned char * _e_hat); 433 434 // decode symbol, returning 0/1/2 for zero/one/multiple errors detected 435 // _sym_enc : input symbol [size: 8 x 1] 436 // _sym_dec : input symbol [size: 9 x 1] 437 int fec_secded7264_decode_symbol(unsigned char * _sym_enc, 438 unsigned char * _sym_dec); 439 440 extern unsigned char secded7264_P[64]; 441 extern unsigned char secded7264_syndrome_w1[72]; 442 443 fec fec_secded7264_create(void *_opts); 444 void fec_secded7264_destroy(fec _q); 445 void fec_secded7264_print(fec _q); 446 void fec_secded7264_encode(fec _q, 447 unsigned int _dec_msg_len, 448 unsigned char * _msg_dec, 449 unsigned char * _msg_enc); 450 void fec_secded7264_decode(fec _q, 451 unsigned int _dec_msg_len, 452 unsigned char * _msg_enc, 453 unsigned char * _msg_dec); 454 455 456 // Convolutional: r1/2 K=7 457 // r1/2 K=9 458 // r1/3 K=9 459 // r1/6 K=15 460 461 // compute encoded message length for block codes 462 // _dec_msg_len : decoded message length (bytes) 463 // _m : input block size (bits) 464 // _k : output block size (bits) 465 unsigned int fec_block_get_enc_msg_len(unsigned int _dec_msg_len, 466 unsigned int _m, 467 unsigned int _k); 468 469 // compute encoded message length for convolutional codes 470 // _dec_msg_len : decoded message length 471 // _K : constraint length 472 // _p : puncturing rate, r = _p / (_p+1) 473 unsigned int fec_conv_get_enc_msg_len(unsigned int _dec_msg_len, 474 unsigned int _K, 475 unsigned int _p); 476 477 // convolutional code polynomials 478 extern int fec_conv27_poly[2]; 479 extern int fec_conv29_poly[2]; 480 extern int fec_conv39_poly[3]; 481 extern int fec_conv615_poly[6]; 482 483 // convolutional code puncturing matrices [R x P] 484 extern int fec_conv27p23_matrix[4]; // [2 x 2] 485 extern int fec_conv27p34_matrix[6]; // [2 x 3] 486 extern int fec_conv27p45_matrix[8]; // [2 x 4] 487 extern int fec_conv27p56_matrix[10]; // [2 x 5] 488 extern int fec_conv27p67_matrix[12]; // [2 x 6] 489 extern int fec_conv27p78_matrix[14]; // [2 x 7] 490 491 extern int fec_conv29p23_matrix[4]; // [2 x 2] 492 extern int fec_conv29p34_matrix[6]; // [2 x 3] 493 extern int fec_conv29p45_matrix[8]; // [2 x 4] 494 extern int fec_conv29p56_matrix[10]; // [2 x 5] 495 extern int fec_conv29p67_matrix[12]; // [2 x 6] 496 extern int fec_conv29p78_matrix[14]; // [2 x 7] 497 498 fec fec_conv_create(fec_scheme _fs); 499 void fec_conv_destroy(fec _q); 500 void fec_conv_print(fec _q); 501 void fec_conv_encode(fec _q, 502 unsigned int _dec_msg_len, 503 unsigned char * _msg_dec, 504 unsigned char * _msg_enc); 505 void fec_conv_decode_hard(fec _q, 506 unsigned int _dec_msg_len, 507 unsigned char * _msg_enc, 508 unsigned char * _msg_dec); 509 void fec_conv_decode_soft(fec _q, 510 unsigned int _dec_msg_len, 511 unsigned char * _msg_enc, 512 unsigned char * _msg_dec); 513 void fec_conv_decode(fec _q, 514 unsigned char * _msg_dec); 515 void fec_conv_setlength(fec _q, 516 unsigned int _dec_msg_len); 517 518 // internal initialization methods (sets r, K, viterbi methods) 519 void fec_conv_init_v27(fec _q); 520 void fec_conv_init_v29(fec _q); 521 void fec_conv_init_v39(fec _q); 522 void fec_conv_init_v615(fec _q); 523 524 // punctured convolutional codes 525 fec fec_conv_punctured_create(fec_scheme _fs); 526 void fec_conv_punctured_destroy(fec _q); 527 void fec_conv_punctured_print(fec _q); 528 void fec_conv_punctured_encode(fec _q, 529 unsigned int _dec_msg_len, 530 unsigned char * _msg_dec, 531 unsigned char * _msg_enc); 532 void fec_conv_punctured_decode_hard(fec _q, 533 unsigned int _dec_msg_len, 534 unsigned char * _msg_enc, 535 unsigned char * _msg_dec); 536 void fec_conv_punctured_decode_soft(fec _q, 537 unsigned int _dec_msg_len, 538 unsigned char * _msg_enc, 539 unsigned char * _msg_dec); 540 void fec_conv_punctured_setlength(fec _q, 541 unsigned int _dec_msg_len); 542 543 // internal initialization methods (sets r, K, viterbi methods, 544 // and puncturing matrix) 545 void fec_conv_init_v27p23(fec _q); 546 void fec_conv_init_v27p34(fec _q); 547 void fec_conv_init_v27p45(fec _q); 548 void fec_conv_init_v27p56(fec _q); 549 void fec_conv_init_v27p67(fec _q); 550 void fec_conv_init_v27p78(fec _q); 551 552 void fec_conv_init_v29p23(fec _q); 553 void fec_conv_init_v29p34(fec _q); 554 void fec_conv_init_v29p45(fec _q); 555 void fec_conv_init_v29p56(fec _q); 556 void fec_conv_init_v29p67(fec _q); 557 void fec_conv_init_v29p78(fec _q); 558 559 // Reed-Solomon 560 561 // compute encoded message length for Reed-Solomon codes 562 // _dec_msg_len : decoded message length 563 // _nroots : number of roots in polynomial 564 // _nn : 565 // _kk : 566 unsigned int fec_rs_get_enc_msg_len(unsigned int _dec_msg_len, 567 unsigned int _nroots, 568 unsigned int _nn, 569 unsigned int _kk); 570 571 572 fec fec_rs_create(fec_scheme _fs); 573 void fec_rs_destroy(fec _q); 574 void fec_rs_init_p8(fec _q); 575 void fec_rs_setlength(fec _q, 576 unsigned int _dec_msg_len); 577 void fec_rs_encode(fec _q, 578 unsigned int _dec_msg_len, 579 unsigned char * _msg_dec, 580 unsigned char * _msg_enc); 581 void fec_rs_decode(fec _q, 582 unsigned int _dec_msg_len, 583 unsigned char * _msg_enc, 584 unsigned char * _msg_dec); 585 586 // phi(x) = -logf( tanhf( x/2 ) ) 587 float sumproduct_phi(float _x); 588 589 // iterate over the sum-product algorithm: 590 // returns 1 if parity checks, 0 otherwise 591 // _m : rows 592 // _n : cols 593 // _H : sparse binary parity check matrix [size: _m x _n] 594 // _LLR : received signal (soft bits, LLR) [size: _n x 1] 595 // _c_hat : estimated transmitted signal [size: _n x 1] 596 // _max_steps : maximum number of steps before bailing 597 int fec_sumproduct(unsigned int _m, 598 unsigned int _n, 599 smatrixb _H, 600 float * _LLR, 601 unsigned char * _c_hat, 602 unsigned int _max_steps); 603 604 // sum-product algorithm, returns 1 if parity checks, 0 otherwise 605 // _m : rows 606 // _n : cols 607 // _H : sparse binary parity check matrix [size: _m x _n] 608 // _c_hat : estimated transmitted signal [size: _n x 1] 609 // 610 // internal state arrays 611 // _Lq : [size: _m x _n] 612 // _Lr : [size: _m x _n] 613 // _Lc : [size: _n x 1] 614 // _LQ : [size: _n x 1] 615 // _parity : _H * _c_hat [size: _m x 1] 616 int fec_sumproduct_step(unsigned int _m, 617 unsigned int _n, 618 smatrixb _H, 619 unsigned char * _c_hat, 620 float * _Lq, 621 float * _Lr, 622 float * _Lc, 623 float * _LQ, 624 unsigned char * _parity); 625 626 // 627 // packetizer 628 // 629 630 // fec/interleaver plan 631 struct fecintlv_plan { 632 unsigned int dec_msg_len; 633 unsigned int enc_msg_len; 634 635 // fec codec 636 fec_scheme fs; 637 fec f; 638 639 // interleaver 640 interleaver q; 641 }; 642 643 #define PACKETIZER_VERSION (1) 644 645 // packetizer object 646 struct packetizer_s { 647 unsigned int msg_len; 648 unsigned int packet_len; 649 650 crc_scheme check; 651 unsigned int crc_length; 652 653 struct fecintlv_plan * plan; 654 unsigned int plan_len; 655 656 // buffers (ping-pong) 657 unsigned int buffer_len; 658 unsigned char * buffer_0; 659 unsigned char * buffer_1; 660 }; 661 662 663 // 664 // MODULE : fft (fast discrete Fourier transform) 665 // 666 667 // fast fourier transform method 668 typedef enum { 669 LIQUID_FFT_METHOD_UNKNOWN=0, // unknown method 670 LIQUID_FFT_METHOD_RADIX2, // Radix-2 (decimation in time) 671 LIQUID_FFT_METHOD_MIXED_RADIX, // Cooley-Tukey mixed-radix FFT (decimation in time) 672 LIQUID_FFT_METHOD_RADER, // Rader's method for FFTs of prime length 673 LIQUID_FFT_METHOD_RADER2, // Rader's method for FFTs of prime length (alternate) 674 LIQUID_FFT_METHOD_DFT, // regular discrete Fourier transform 675 } liquid_fft_method; 676 677 // Macro : FFT (internal) 678 // FFT : name-mangling macro 679 // T : primitive data type 680 // TC : primitive data type (complex) 681 #define LIQUID_FFT_DEFINE_INTERNAL_API(FFT,T,TC) \ 682 \ 683 /* print plan recursively */ \ 684 void FFT(_print_plan_recursive)(FFT(plan) _q, \ 685 unsigned int _level); \ 686 \ 687 /* type definitions for create/destroy/execute functions */ \ 688 typedef FFT(plan)(FFT(_create_t)) (unsigned int _nfft, \ 689 TC * _x, \ 690 TC * _y, \ 691 int _dir, \ 692 int _flags); \ 693 typedef void (FFT(_destroy_t))(FFT(plan) _q); \ 694 typedef void (FFT(_execute_t))(FFT(plan) _q); \ 695 \ 696 /* FFT create methods */ \ 697 FFT(_create_t) FFT(_create_plan_dft); \ 698 FFT(_create_t) FFT(_create_plan_radix2); \ 699 FFT(_create_t) FFT(_create_plan_mixed_radix); \ 700 FFT(_create_t) FFT(_create_plan_rader); \ 701 FFT(_create_t) FFT(_create_plan_rader2); \ 702 \ 703 /* FFT destroy methods */ \ 704 FFT(_destroy_t) FFT(_destroy_plan_dft); \ 705 FFT(_destroy_t) FFT(_destroy_plan_radix2); \ 706 FFT(_destroy_t) FFT(_destroy_plan_mixed_radix); \ 707 FFT(_destroy_t) FFT(_destroy_plan_rader); \ 708 FFT(_destroy_t) FFT(_destroy_plan_rader2); \ 709 \ 710 /* FFT execute methods */ \ 711 FFT(_execute_t) FFT(_execute_dft); \ 712 FFT(_execute_t) FFT(_execute_radix2); \ 713 FFT(_execute_t) FFT(_execute_mixed_radix); \ 714 FFT(_execute_t) FFT(_execute_rader); \ 715 FFT(_execute_t) FFT(_execute_rader2); \ 716 \ 717 /* specific codelets for small DFTs */ \ 718 FFT(_execute_t) FFT(_execute_dft_2); \ 719 FFT(_execute_t) FFT(_execute_dft_3); \ 720 FFT(_execute_t) FFT(_execute_dft_4); \ 721 FFT(_execute_t) FFT(_execute_dft_5); \ 722 FFT(_execute_t) FFT(_execute_dft_6); \ 723 FFT(_execute_t) FFT(_execute_dft_7); \ 724 FFT(_execute_t) FFT(_execute_dft_8); \ 725 FFT(_execute_t) FFT(_execute_dft_16); \ 726 \ 727 /* additional methods */ \ 728 unsigned int FFT(_estimate_mixed_radix)(unsigned int _nfft); \ 729 \ 730 /* discrete cosine transform (DCT) prototypes */ \ 731 void FFT(_execute_REDFT00)(FFT(plan) _q); /* DCT-I */ \ 732 void FFT(_execute_REDFT10)(FFT(plan) _q); /* DCT-II */ \ 733 void FFT(_execute_REDFT01)(FFT(plan) _q); /* DCT-III */ \ 734 void FFT(_execute_REDFT11)(FFT(plan) _q); /* DCT-IV */ \ 735 \ 736 /* discrete sine transform (DST) prototypes */ \ 737 void FFT(_execute_RODFT00)(FFT(plan) _q); /* DST-I */ \ 738 void FFT(_execute_RODFT10)(FFT(plan) _q); /* DST-II */ \ 739 void FFT(_execute_RODFT01)(FFT(plan) _q); /* DST-III */ \ 740 void FFT(_execute_RODFT11)(FFT(plan) _q); /* DST-IV */ \ 741 \ 742 /* destroy real-to-real one-dimensional plan */ \ 743 void FFT(_destroy_plan_r2r_1d)(FFT(plan) _q); \ 744 \ 745 /* print real-to-real one-dimensional plan */ \ 746 void FFT(_print_plan_r2r_1d)(FFT(plan) _q); \ 747 748 // determine best FFT method based on size 749 liquid_fft_method liquid_fft_estimate_method(unsigned int _nfft); 750 751 // is input radix-2? 752 int fft_is_radix2(unsigned int _n); 753 754 // miscellaneous functions 755 unsigned int fft_reverse_index(unsigned int _i, unsigned int _n); 756 757 758 LIQUID_FFT_DEFINE_INTERNAL_API(LIQUID_FFT_MANGLE_FLOAT, float, liquid_float_complex) 759 760 // Use fftw library if installed (and not overridden with configuration), 761 // otherwise use internal (less efficient) fft library. 762 #if HAVE_FFTW3_H && !defined LIQUID_FFTOVERRIDE 763 # include <fftw3.h> 764 # define FFT_PLAN fftwf_plan 765 # define FFT_CREATE_PLAN fftwf_plan_dft_1d 766 # define FFT_DESTROY_PLAN fftwf_destroy_plan 767 # define FFT_EXECUTE fftwf_execute 768 # define FFT_DIR_FORWARD FFTW_FORWARD 769 # define FFT_DIR_BACKWARD FFTW_BACKWARD 770 # define FFT_METHOD FFTW_ESTIMATE 771 #else 772 # define FFT_PLAN fftplan 773 # define FFT_CREATE_PLAN fft_create_plan 774 # define FFT_DESTROY_PLAN fft_destroy_plan 775 # define FFT_EXECUTE fft_execute 776 # define FFT_DIR_FORWARD LIQUID_FFT_FORWARD 777 # define FFT_DIR_BACKWARD LIQUID_FFT_BACKWARD 778 # define FFT_METHOD 0 779 #endif 780 781 782 783 // 784 // MODULE : filter 785 // 786 787 // esimate required filter length given transition bandwidth and 788 // stop-band attenuation (algorithm from [Vaidyanathan:1993]) 789 // _df : transition bandwidth (0 < _df < 0.5) 790 // _As : stop-band attenuation [dB] (As > 0) 791 float estimate_req_filter_len_Kaiser(float _df, 792 float _As); 793 794 // esimate required filter length given transition bandwidth and 795 // stop-band attenuation (algorithm from [Herrmann:1973]) 796 // _df : transition bandwidth (0 < _df < 0.5) 797 // _As : stop-band attenuation [dB] (As > 0) 798 float estimate_req_filter_len_Herrmann(float _df, 799 float _As); 800 801 802 // fir_farrow 803 #define LIQUID_FIRFARROW_DEFINE_INTERNAL_API(FIRFARROW,TO,TC,TI) \ 804 void FIRFARROW(_genpoly)(FIRFARROW() _q); 805 806 LIQUID_FIRFARROW_DEFINE_INTERNAL_API(LIQUID_FIRFARROW_MANGLE_RRRF, 807 float, 808 float, 809 float) 810 811 LIQUID_FIRFARROW_DEFINE_INTERNAL_API(LIQUID_FIRFARROW_MANGLE_CRCF, 812 liquid_float_complex, 813 float, 814 liquid_float_complex) 815 816 817 818 // 819 // iirfiltsos : infinite impulse respone filter (second-order sections) 820 // 821 #define LIQUID_IIRFILTSOS_MANGLE_RRRF(name) LIQUID_CONCAT(iirfiltsos_rrrf,name) 822 #define LIQUID_IIRFILTSOS_MANGLE_CRCF(name) LIQUID_CONCAT(iirfiltsos_crcf,name) 823 #define LIQUID_IIRFILTSOS_MANGLE_CCCF(name) LIQUID_CONCAT(iirfiltsos_cccf,name) 824 825 #define LIQUID_IIRFILTSOS_DEFINE_INTERNAL_API(IIRFILTSOS,TO,TC,TI) \ 826 typedef struct IIRFILTSOS(_s) * IIRFILTSOS(); \ 827 \ 828 /* create 2nd-order infinite impulse reponse filter */ \ 829 /* _b : feed-forward coefficients [size: _3 x 1] */ \ 830 /* _a : feed-back coefficients [size: _3 x 1] */ \ 831 IIRFILTSOS() IIRFILTSOS(_create)(TC * _b, \ 832 TC * _a); \ 833 \ 834 /* explicitly set 2nd-order IIR filter coefficients */ \ 835 /* _q : iirfiltsos object */ \ 836 /* _b : feed-forward coefficients [size: _3 x 1] */ \ 837 /* _a : feed-back coefficients [size: _3 x 1] */ \ 838 void IIRFILTSOS(_set_coefficients)(IIRFILTSOS() _q, \ 839 TC * _b, \ 840 TC * _a); \ 841 \ 842 /* destroy iirfiltsos object, freeing all internal memory */ \ 843 void IIRFILTSOS(_destroy)(IIRFILTSOS() _q); \ 844 \ 845 /* print iirfiltsos object properties to stdout */ \ 846 void IIRFILTSOS(_print)(IIRFILTSOS() _q); \ 847 \ 848 /* clear/reset iirfiltsos object internals */ \ 849 void IIRFILTSOS(_reset)(IIRFILTSOS() _q); \ 850 \ 851 /* compute filter output */ \ 852 /* _q : iirfiltsos object */ \ 853 /* _x : input sample */ \ 854 /* _y : output sample pointer */ \ 855 void IIRFILTSOS(_execute)(IIRFILTSOS() _q, \ 856 TI _x, \ 857 TO * _y); \ 858 \ 859 /* compute filter output, direct-form I method */ \ 860 /* _q : iirfiltsos object */ \ 861 /* _x : input sample */ \ 862 /* _y : output sample pointer */ \ 863 void IIRFILTSOS(_execute_df1)(IIRFILTSOS() _q, \ 864 TI _x, \ 865 TO * _y); \ 866 \ 867 /* compute filter output, direct-form II method */ \ 868 /* _q : iirfiltsos object */ \ 869 /* _x : input sample */ \ 870 /* _y : output sample pointer */ \ 871 void IIRFILTSOS(_execute_df2)(IIRFILTSOS() _q, \ 872 TI _x, \ 873 TO * _y); \ 874 \ 875 /* compute and return group delay of filter object */ \ 876 /* _q : filter object */ \ 877 /* _fc : frequency to evaluate */ \ 878 float IIRFILTSOS(_groupdelay)(IIRFILTSOS() _q, \ 879 float _fc); \ 880 881 LIQUID_IIRFILTSOS_DEFINE_INTERNAL_API(LIQUID_IIRFILTSOS_MANGLE_RRRF, 882 float, 883 float, 884 float) 885 886 LIQUID_IIRFILTSOS_DEFINE_INTERNAL_API(LIQUID_IIRFILTSOS_MANGLE_CRCF, 887 liquid_float_complex, 888 float, 889 liquid_float_complex) 890 891 LIQUID_IIRFILTSOS_DEFINE_INTERNAL_API(LIQUID_IIRFILTSOS_MANGLE_CCCF, 892 liquid_float_complex, 893 liquid_float_complex, 894 liquid_float_complex) 895 896 897 // firdes : finite impulse response filter design 898 899 // Find approximate bandwidth adjustment factor rho based on 900 // filter delay and desired excess bandwdith factor. 901 // 902 // _m : filter delay (symbols) 903 // _beta : filter excess bandwidth factor (0,1) 904 float rkaiser_approximate_rho(unsigned int _m, 905 float _beta); 906 907 // Design frequency-shifted root-Nyquist filter based on 908 // the Kaiser-windowed sinc using the bisection method 909 // 910 // _k : filter over-sampling rate (samples/symbol) 911 // _m : filter delay (symbols) 912 // _beta : filter excess bandwidth factor (0,1) 913 // _dt : filter fractional sample delay 914 // _h : resulting filter [size: 2*_k*_m+1] 915 // _rho : transition bandwidth adjustment, 0 < _rho < 1 916 void liquid_firdes_rkaiser_bisection(unsigned int _k, 917 unsigned int _m, 918 float _beta, 919 float _dt, 920 float * _h, 921 float * _rho); 922 923 // Design frequency-shifted root-Nyquist filter based on 924 // the Kaiser-windowed sinc using the quadratic method. 925 // 926 // _k : filter over-sampling rate (samples/symbol) 927 // _m : filter delay (symbols) 928 // _beta : filter excess bandwidth factor (0,1) 929 // _dt : filter fractional sample delay 930 // _h : resulting filter [size: 2*_k*_m+1] 931 // _rho : transition bandwidth adjustment, 0 < _rho < 1 932 void liquid_firdes_rkaiser_quadratic(unsigned int _k, 933 unsigned int _m, 934 float _beta, 935 float _dt, 936 float * _h, 937 float * _rho); 938 939 // compute filter coefficients and determine resulting ISI 940 // 941 // _k : filter over-sampling rate (samples/symbol) 942 // _m : filter delay (symbols) 943 // _beta : filter excess bandwidth factor (0,1) 944 // _dt : filter fractional sample delay 945 // _rho : transition bandwidth adjustment, 0 < _rho < 1 946 // _h : filter buffer [size: 2*_k*_m+1] 947 float liquid_firdes_rkaiser_internal_isi(unsigned int _k, 948 unsigned int _m, 949 float _beta, 950 float _dt, 951 float _rho, 952 float * _h); 953 954 // Design flipped Nyquist/root-Nyquist filters 955 void liquid_firdes_fnyquist(liquid_firfilt_type _type, 956 int _root, 957 unsigned int _k, 958 unsigned int _m, 959 float _beta, 960 float _dt, 961 float * _h); 962 963 // flipped exponential frequency response 964 void liquid_firdes_fexp_freqresponse(unsigned int _k, 965 unsigned int _m, 966 float _beta, 967 float * _H); 968 969 // flipped hyperbolic secant frequency response 970 void liquid_firdes_fsech_freqresponse(unsigned int _k, 971 unsigned int _m, 972 float _beta, 973 float * _H); 974 975 // flipped hyperbolic secant frequency response 976 void liquid_firdes_farcsech_freqresponse(unsigned int _k, 977 unsigned int _m, 978 float _beta, 979 float * _H); 980 981 // iirdes : infinite impulse response filter design 982 983 // Sorts array _z of complex numbers into complex conjugate pairs to 984 // within a tolerance. Conjugate pairs are ordered by increasing real 985 // component with the negative imaginary element first. All pure-real 986 // elements are placed at the end of the array. 987 // 988 // Example: 989 // v: liquid_cplxpair(v): 990 // 10 + j*3 -3 - j*4 991 // 5 + j*0 3 + j*4 992 // -3 + j*4 10 - j*3 993 // 10 - j*3 10 + j*3 994 // 3 + j*0 3 + j*0 995 // -3 + j*4 5 + j*0 996 // 997 // _z : complex array (size _n) 998 // _n : number of elements in _z 999 // _tol : tolerance for finding complex pairs 1000 // _p : resulting pairs, pure real values of _z at end 1001 void liquid_cplxpair(float complex * _z, 1002 unsigned int _n, 1003 float _tol, 1004 float complex * _p); 1005 1006 // post-process cleanup used with liquid_cplxpair 1007 // 1008 // once pairs have been identified and ordered, this method 1009 // will clean up the result by ensuring the following: 1010 // * pairs are perfect conjugates 1011 // * pairs have negative imaginary component first 1012 // * pairs are ordered by increasing real component 1013 // * pure-real elements are ordered by increasing value 1014 // 1015 // _p : pre-processed complex array [size: _n x 1] 1016 // _n : array length 1017 // _num_pairs : number of complex conjugate pairs 1018 void liquid_cplxpair_cleanup(float complex * _p, 1019 unsigned int _n, 1020 unsigned int _num_pairs); 1021 1022 // Jacobian elliptic functions (src/filter/src/ellip.c) 1023 1024 // Landen transformation (_n iterations) 1025 void landenf(float _k, 1026 unsigned int _n, 1027 float * _v); 1028 1029 // compute elliptic integral K(k) for _n recursions 1030 void ellipkf(float _k, 1031 unsigned int _n, 1032 float * _K, 1033 float * _Kp); 1034 1035 // elliptic degree 1036 float ellipdegf(float _N, 1037 float _k1, 1038 unsigned int _n); 1039 1040 // elliptic cd() function (_n recursions) 1041 float complex ellip_cdf(float complex _u, 1042 float _k, 1043 unsigned int _n); 1044 1045 // elliptic inverse cd() function (_n recursions) 1046 float complex ellip_acdf(float complex _u, 1047 float _k, 1048 unsigned int _n); 1049 1050 // elliptic sn() function (_n recursions) 1051 float complex ellip_snf(float complex _u, 1052 float _k, 1053 unsigned int _n); 1054 1055 // elliptic inverse sn() function (_n recursions) 1056 float complex ellip_asnf(float complex _u, 1057 float _k, 1058 unsigned int _n); 1059 1060 // 1061 // MODULE : framing 1062 // 1063 1064 // 1065 // bpacket 1066 // 1067 1068 #define BPACKET_VERSION (101+PACKETIZER_VERSION) 1069 1070 // generator 1071 void bpacketgen_compute_packet_len(bpacketgen _q); 1072 void bpacketgen_assemble_pnsequence(bpacketgen _q); 1073 void bpacketgen_assemble_header(bpacketgen _q); 1074 1075 // synchronizer 1076 void bpacketsync_assemble_pnsequence(bpacketsync _q); 1077 void bpacketsync_execute_seekpn(bpacketsync _q, unsigned char _bit); 1078 void bpacketsync_execute_rxheader(bpacketsync _q, unsigned char _bit); 1079 void bpacketsync_execute_rxpayload(bpacketsync _q, unsigned char _bit); 1080 void bpacketsync_decode_header(bpacketsync _q); 1081 void bpacketsync_decode_payload(bpacketsync _q); 1082 void bpacketsync_reconfig(bpacketsync _q); 1083 1084 1085 // 1086 // flexframe 1087 // 1088 1089 // flexframe protocol 1090 #define FLEXFRAME_PROTOCOL (101+PACKETIZER_VERSION) 1091 1092 // header description 1093 #define FLEXFRAME_H_USER_DEFAULT (14) // default length for user-defined array 1094 #define FLEXFRAME_H_DEC (6) // decoded length 1095 #define FLEXFRAME_H_CRC (LIQUID_CRC_32) // header CRC 1096 #define FLEXFRAME_H_FEC0 (LIQUID_FEC_SECDED7264) // header FEC (inner) 1097 #define FLEXFRAME_H_FEC1 (LIQUID_FEC_HAMMING84) // header FEC (outer) 1098 #define FLEXFRAME_H_MOD (LIQUID_MODEM_QPSK) // modulation scheme 1099 1100 1101 // 1102 // gmskframe 1103 // 1104 1105 #define GMSKFRAME_VERSION (3+PACKETIZER_VERSION) 1106 1107 // header description 1108 #define GMSKFRAME_H_USER_DEFAULT (8) // user-defined array 1109 #define GMSKFRAME_H_DEC (5) // decoded length 1110 #define GMSKFRAME_H_CRC (LIQUID_CRC_32) // header CRC 1111 #define GMSKFRAME_H_FEC (LIQUID_FEC_HAMMING128) // header FEC 1112 1113 1114 // 1115 // ofdmflexframe 1116 // 1117 1118 #define OFDMFLEXFRAME_PROTOCOL (104+PACKETIZER_VERSION) 1119 1120 // header description 1121 #define OFDMFLEXFRAME_H_USER_DEFAULT (8) // default length for user-defined array 1122 #define OFDMFLEXFRAME_H_DEC (6) // decoded length 1123 #define OFDMFLEXFRAME_H_CRC (LIQUID_CRC_32) // header CRC 1124 #define OFDMFLEXFRAME_H_FEC0 (LIQUID_FEC_GOLAY2412) // header FEC (inner) 1125 #define OFDMFLEXFRAME_H_FEC1 (LIQUID_FEC_NONE) // header FEC (outer) 1126 #define OFDMFLEXFRAME_H_MOD (LIQUID_MODEM_BPSK) // modulation scheme 1127 1128 1129 // 1130 // dsssframe 1131 // 1132 1133 #define DSSSFRAME_PROTOCOL (101 + PACKETIZER_VERSION) 1134 #define DSSSFRAME_H_USER_DEFAULT (8) 1135 #define DSSSFRAME_H_DEC (5) 1136 #define DSSSFRAME_H_CRC (LIQUID_CRC_32) 1137 #define DSSSFRAME_H_FEC0 (LIQUID_FEC_GOLAY2412) 1138 #define DSSSFRAME_H_FEC1 (LIQUID_FEC_NONE) 1139 1140 // 1141 // multi-signal source for testing (no meaningful data, just signals) 1142 // 1143 #define LIQUID_QSOURCE_MANGLE_CFLOAT(name) LIQUID_CONCAT(qsourcecf,name) 1144 1145 #define LIQUID_QSOURCE_DEFINE_API(QSOURCE,TO) \ 1146 \ 1147 /* Multi-signal source generator object */ \ 1148 typedef struct QSOURCE(_s) * QSOURCE(); \ 1149 \ 1150 /* Create default qsource object, type uninitialized */ \ 1151 QSOURCE() QSOURCE(_create)(unsigned int _M, \ 1152 unsigned int _m, \ 1153 float _As, \ 1154 float _fc, \ 1155 float _bw, \ 1156 float _gain); \ 1157 \ 1158 /* Initialize user-defined qsource object */ \ 1159 void QSOURCE(_init_user)(QSOURCE() _q, \ 1160 void * _userdata, \ 1161 void * _callback); \ 1162 \ 1163 /* Initialize qsource tone object */ \ 1164 void QSOURCE(_init_tone)(QSOURCE() _q); \ 1165 \ 1166 /* Add chirp to signal generator, returning id of signal */ \ 1167 /* _q : signal source object */ \ 1168 /* _duration : duration of chirp [samples] */ \ 1169 /* _negate : negate frequency direction */ \ 1170 /* _repeat : repeat signal? or just run once */ \ 1171 void QSOURCE(_init_chirp)(QSOURCE() _q, \ 1172 float _duration, \ 1173 int _negate, \ 1174 int _repeat); \ 1175 \ 1176 /* Initialize qsource noise object */ \ 1177 void QSOURCE(_init_noise)(QSOURCE() _q); \ 1178 \ 1179 /* Initialize qsource linear modem object */ \ 1180 void QSOURCE(_init_modem)(QSOURCE() _q, \ 1181 int _ms, \ 1182 unsigned int _m, \ 1183 float _beta); \ 1184 \ 1185 /* Initialize frequency-shift keying modem signal source */ \ 1186 /* _q : signal source object */ \ 1187 /* _m : bits per symbol, _bps > 0 */ \ 1188 /* _k : samples/symbol, _k >= 2^_m */ \ 1189 void QSOURCE(_init_fsk)(QSOURCE() _q, \ 1190 unsigned int _m, \ 1191 unsigned int _k); \ 1192 \ 1193 /* Initialize qsource GMSK modem object */ \ 1194 void QSOURCE(_init_gmsk)(QSOURCE() _q, \ 1195 unsigned int _m, \ 1196 float _bt); \ 1197 \ 1198 /* Destroy qsource object */ \ 1199 void QSOURCE(_destroy)(QSOURCE() _q); \ 1200 \ 1201 /* Print qsource object */ \ 1202 void QSOURCE(_print)(QSOURCE() _q); \ 1203 \ 1204 /* Reset qsource object */ \ 1205 void QSOURCE(_reset)(QSOURCE() _q); \ 1206 \ 1207 /* Get/set source id */ \ 1208 void QSOURCE(_set_id)(QSOURCE() _q, int _id); \ 1209 int QSOURCE(_get_id)(QSOURCE() _q); \ 1210 \ 1211 void QSOURCE(_enable)(QSOURCE() _q); \ 1212 void QSOURCE(_disable)(QSOURCE() _q); \ 1213 \ 1214 void QSOURCE(_set_gain)(QSOURCE() _q, \ 1215 float _gain_dB); \ 1216 \ 1217 float QSOURCE(_get_gain)(QSOURCE() _q); \ 1218 \ 1219 /* Get number of samples generated by the object so far */ \ 1220 /* _q : msource object */ \ 1221 /* _gain : signal gain output [dB] */ \ 1222 uint64_t QSOURCE(_get_num_samples)(QSOURCE() _q); \ 1223 \ 1224 void QSOURCE(_set_frequency)(QSOURCE() _q, \ 1225 float _dphi); \ 1226 \ 1227 float QSOURCE(_get_frequency)(QSOURCE() _q); \ 1228 \ 1229 void QSOURCE(_generate)(QSOURCE() _q, \ 1230 TO * _v); \ 1231 \ 1232 void QSOURCE(_generate_into)(QSOURCE() _q, \ 1233 TO * _buf); \ 1234 1235 LIQUID_QSOURCE_DEFINE_API(LIQUID_QSOURCE_MANGLE_CFLOAT, liquid_float_complex) 1236 1237 // 1238 // MODULE : math 1239 // 1240 1241 // 1242 // basic trigonometric functions 1243 // 1244 float liquid_sinf(float _x); 1245 float liquid_cosf(float _x); 1246 float liquid_tanf(float _x); 1247 void liquid_sincosf(float _x, 1248 float * _sinf, 1249 float * _cosf); 1250 float liquid_expf(float _x); 1251 float liquid_logf(float _x); 1252 1253 // 1254 // complex math operations 1255 // 1256 1257 // complex square root 1258 float complex liquid_csqrtf(float complex _z); 1259 1260 // complex exponent, logarithm 1261 float complex liquid_cexpf(float complex _z); 1262 float complex liquid_clogf(float complex _z); 1263 1264 // complex arcsin, arccos, arctan 1265 float complex liquid_casinf(float complex _z); 1266 float complex liquid_cacosf(float complex _z); 1267 float complex liquid_catanf(float complex _z); 1268 1269 // faster approximation to arg{*} 1270 float liquid_cargf_approx(float complex _z); 1271 1272 1273 // internal trig helper functions 1274 1275 // complex rotation vector: cexpf(_Complex_I*THETA) 1276 #define liquid_cexpjf(THETA) (cosf(THETA) + _Complex_I*sinf(THETA)) 1277 1278 1279 // 1280 // MODULE : matrix 1281 // 1282 1283 // large macro 1284 // MATRIX : name-mangling macro 1285 // T : data type 1286 #define LIQUID_MATRIX_DEFINE_INTERNAL_API(MATRIX,T) \ 1287 T MATRIX(_det2x2)(T * _x, \ 1288 unsigned int _rx, \ 1289 unsigned int _cx); 1290 1291 1292 LIQUID_MATRIX_DEFINE_INTERNAL_API(LIQUID_MATRIX_MANGLE_FLOAT, float) 1293 LIQUID_MATRIX_DEFINE_INTERNAL_API(LIQUID_MATRIX_MANGLE_DOUBLE, double) 1294 1295 LIQUID_MATRIX_DEFINE_INTERNAL_API(LIQUID_MATRIX_MANGLE_CFLOAT, liquid_float_complex) 1296 LIQUID_MATRIX_DEFINE_INTERNAL_API(LIQUID_MATRIX_MANGLE_CDOUBLE, liquid_double_complex) 1297 1298 1299 // sparse 'alist' matrix type (similar to MacKay, Davey Lafferty convention) 1300 // large macro 1301 // SMATRIX : name-mangling macro 1302 // T : primitive data type 1303 #define LIQUID_SMATRIX_DEFINE_INTERNAL_API(SMATRIX,T) \ 1304 \ 1305 void SMATRIX(_reset_max_mlist)(SMATRIX() _q); \ 1306 void SMATRIX(_reset_max_nlist)(SMATRIX() _q); \ 1307 1308 LIQUID_SMATRIX_DEFINE_INTERNAL_API(LIQUID_SMATRIX_MANGLE_BOOL, unsigned char) 1309 LIQUID_SMATRIX_DEFINE_INTERNAL_API(LIQUID_SMATRIX_MANGLE_FLOAT, float) 1310 LIQUID_SMATRIX_DEFINE_INTERNAL_API(LIQUID_SMATRIX_MANGLE_INT, short int) 1311 1312 // search for index placement in list 1313 unsigned short int smatrix_indexsearch(unsigned short int * _list, 1314 unsigned int _num_elements, 1315 unsigned short int _value); 1316 1317 1318 1319 1320 // 1321 // MODULE : modem 1322 // 1323 1324 // 'Square' QAM 1325 #define QAM4_ALPHA (1./sqrt(2)) 1326 #define QAM8_ALPHA (1./sqrt(6)) 1327 #define QAM16_ALPHA (1./sqrt(10)) 1328 #define QAM32_ALPHA (1./sqrt(20)) 1329 #define QAM64_ALPHA (1./sqrt(42)) 1330 #define QAM128_ALPHA (1./sqrt(82)) 1331 #define QAM256_ALPHA (1./sqrt(170)) 1332 #define QAM1024_ALPHA (1./sqrt(682)) 1333 #define QAM4096_ALPHA (1./sqrt(2730)) 1334 1335 // Rectangular QAM 1336 #define RQAM4_ALPHA QAM4_ALPHA 1337 #define RQAM8_ALPHA QAM8_ALPHA 1338 #define RQAM16_ALPHA QAM16_ALPHA 1339 #define RQAM32_ALPHA (1./sqrt(26)) 1340 #define RQAM64_ALPHA QAM64_ALPHA 1341 #define RQAM128_ALPHA (1./sqrt(106)) 1342 #define RQAM256_ALPHA QAM256_ALPHA 1343 #define RQAM512_ALPHA (1./sqrt(426)) 1344 #define RQAM1024_ALPHA QAM1024_ALPHA 1345 #define RQAM2048_ALPHA (1./sqrt(1706)) 1346 #define RQAM4096_ALPHA QAM4096_ALPHA 1347 1348 // ASK 1349 #define ASK2_ALPHA (1.) 1350 #define ASK4_ALPHA (1./sqrt(5)) 1351 #define ASK8_ALPHA (1./sqrt(21)) 1352 #define ASK16_ALPHA (1./sqrt(85)) 1353 #define ASK32_ALPHA (1./sqrt(341)) 1354 #define ASK64_ALPHA (1./sqrt(1365)) 1355 #define ASK128_ALPHA (1./sqrt(5461)) 1356 #define ASK256_ALPHA (1./sqrt(21845)) 1357 1358 // Macro : MODEM 1359 // MODEM : name-mangling macro 1360 // T : primitive data type 1361 // TC : primitive data type (complex) 1362 #define LIQUID_MODEM_DEFINE_INTERNAL_API(MODEM,T,TC) \ 1363 \ 1364 /* initialize a generic modem object */ \ 1365 void MODEM(_init)(MODEM() _q, unsigned int _bits_per_symbol); \ 1366 \ 1367 /* initialize symbol map for fast modulation */ \ 1368 void MODEM(_init_map)(MODEM() _q); \ 1369 \ 1370 /* generic modem create routines */ \ 1371 MODEM() MODEM(_create_ask)( unsigned int _bits_per_symbol); \ 1372 MODEM() MODEM(_create_qam)( unsigned int _bits_per_symbol); \ 1373 MODEM() MODEM(_create_psk)( unsigned int _bits_per_symbol); \ 1374 MODEM() MODEM(_create_dpsk)(unsigned int _bits_per_symbol); \ 1375 MODEM() MODEM(_create_apsk)(unsigned int _bits_per_symbol); \ 1376 MODEM() MODEM(_create_arb)( unsigned int _bits_per_symbol); \ 1377 \ 1378 /* Initialize arbitrary modem constellation */ \ 1379 void MODEM(_arb_init)(MODEM() _q, \ 1380 TC * _symbol_map, \ 1381 unsigned int _len); \ 1382 \ 1383 /* Initialize arb modem constellation from external file */ \ 1384 void MODEM(_arb_init_file)(MODEM() _q, char * _filename); \ 1385 \ 1386 /* specific modem create routines */ \ 1387 MODEM() MODEM(_create_bpsk)(void); \ 1388 MODEM() MODEM(_create_qpsk)(void); \ 1389 MODEM() MODEM(_create_ook)(void); \ 1390 MODEM() MODEM(_create_sqam32)(void); \ 1391 MODEM() MODEM(_create_sqam128)(void); \ 1392 MODEM() MODEM(_create_V29)(void); \ 1393 MODEM() MODEM(_create_arb16opt)(void); \ 1394 MODEM() MODEM(_create_arb32opt)(void); \ 1395 MODEM() MODEM(_create_arb64opt)(void); \ 1396 MODEM() MODEM(_create_arb128opt)(void); \ 1397 MODEM() MODEM(_create_arb256opt)(void); \ 1398 MODEM() MODEM(_create_arb64vt)(void); \ 1399 \ 1400 /* Scale arbitrary modem energy to unity */ \ 1401 void MODEM(_arb_scale)(MODEM() _q); \ 1402 \ 1403 /* Balance I/Q */ \ 1404 void MODEM(_arb_balance_iq)(MODEM() _q); \ 1405 \ 1406 /* modulate using symbol map (look-up table) */ \ 1407 void MODEM(_modulate_map)(MODEM() _q, \ 1408 unsigned int _sym_in, \ 1409 TC * _y); \ 1410 \ 1411 /* modem modulate routines */ \ 1412 void MODEM(_modulate_ask) ( MODEM(), unsigned int, TC *); \ 1413 void MODEM(_modulate_qam) ( MODEM(), unsigned int, TC *); \ 1414 void MODEM(_modulate_psk) ( MODEM(), unsigned int, TC *); \ 1415 void MODEM(_modulate_dpsk) ( MODEM(), unsigned int, TC *); \ 1416 void MODEM(_modulate_arb) ( MODEM(), unsigned int, TC *); \ 1417 void MODEM(_modulate_apsk) ( MODEM(), unsigned int, TC *); \ 1418 void MODEM(_modulate_bpsk) ( MODEM(), unsigned int, TC *); \ 1419 void MODEM(_modulate_qpsk) ( MODEM(), unsigned int, TC *); \ 1420 void MODEM(_modulate_ook) ( MODEM(), unsigned int, TC *); \ 1421 void MODEM(_modulate_sqam32) ( MODEM(), unsigned int, TC *); \ 1422 void MODEM(_modulate_sqam128) ( MODEM(), unsigned int, TC *); \ 1423 \ 1424 /* modem demodulate routines */ \ 1425 void MODEM(_demodulate_ask) ( MODEM(), TC, unsigned int *); \ 1426 void MODEM(_demodulate_qam) ( MODEM(), TC, unsigned int *); \ 1427 void MODEM(_demodulate_psk) ( MODEM(), TC, unsigned int *); \ 1428 void MODEM(_demodulate_dpsk) ( MODEM(), TC, unsigned int *); \ 1429 void MODEM(_demodulate_arb) ( MODEM(), TC, unsigned int *); \ 1430 void MODEM(_demodulate_apsk) ( MODEM(), TC, unsigned int *); \ 1431 void MODEM(_demodulate_bpsk) ( MODEM(), TC, unsigned int *); \ 1432 void MODEM(_demodulate_qpsk) ( MODEM(), TC, unsigned int *); \ 1433 void MODEM(_demodulate_ook) ( MODEM(), TC, unsigned int *); \ 1434 void MODEM(_demodulate_sqam32) ( MODEM(), TC, unsigned int *); \ 1435 void MODEM(_demodulate_sqam128)( MODEM(), TC, unsigned int *); \ 1436 \ 1437 /* modem demodulate (soft) routines */ \ 1438 void MODEM(_demodulate_soft_bpsk)(MODEM() _q, \ 1439 TC _x, \ 1440 unsigned int * _sym_out, \ 1441 unsigned char * _soft_bits); \ 1442 void MODEM(_demodulate_soft_qpsk)(MODEM() _q, \ 1443 TC _x, \ 1444 unsigned int * _sym_out, \ 1445 unsigned char * _soft_bits); \ 1446 void MODEM(_demodulate_soft_arb)( MODEM() _q, \ 1447 TC _x, \ 1448 unsigned int * _sym_out, \ 1449 unsigned char * _soft_bits); \ 1450 \ 1451 /* generate soft demodulation look-up table */ \ 1452 void MODEM(_demodsoft_gentab)(MODEM() _q, \ 1453 unsigned int _p); \ 1454 \ 1455 /* generic soft demodulation routine using nearest-neighbors */ \ 1456 /* look-up table */ \ 1457 void MODEM(_demodulate_soft_table)(MODEM() _q, \ 1458 TC _x, \ 1459 unsigned int * _sym_out, \ 1460 unsigned char * _soft_bits); \ 1461 \ 1462 /* Demodulate a linear symbol constellation using dynamic */ \ 1463 /* threshold calculation */ \ 1464 /* _v : input value */ \ 1465 /* _m : bits per symbol */ \ 1466 /* _alpha : scaling factor */ \ 1467 /* _s : demodulated symbol */ \ 1468 /* _res : residual */ \ 1469 void MODEM(_demodulate_linear_array)(T _v, \ 1470 unsigned int _m, \ 1471 T _alpha, \ 1472 unsigned int * _s, \ 1473 T * _res); \ 1474 \ 1475 /* Demodulate a linear symbol constellation using */ \ 1476 /* refereneced lookup table */ \ 1477 /* _v : input value */ \ 1478 /* _m : bits per symbol */ \ 1479 /* _ref : array of thresholds */ \ 1480 /* _s : demodulated symbol */ \ 1481 /* _res : residual */ \ 1482 void MODEM(_demodulate_linear_array_ref)(T _v, \ 1483 unsigned int _m, \ 1484 T * _ref, \ 1485 unsigned int * _s, \ 1486 T * _res); \ 1487 1488 1489 1490 // define internal modem APIs 1491 LIQUID_MODEM_DEFINE_INTERNAL_API(LIQUID_MODEM_MANGLE_FLOAT,float,float complex) 1492 1493 // APSK constants (container for apsk structure definitions) 1494 struct liquid_apsk_s { 1495 modulation_scheme scheme; // APSK modulation scheme 1496 unsigned int num_levels; // number of rings 1497 unsigned int * p; // number of symbols per ring 1498 float * r; // radius of each ring 1499 float * phi; // phase offset of each ring 1500 float * r_slicer; // radius slicer 1501 unsigned char * map; // symbol mapping 1502 }; 1503 1504 extern struct liquid_apsk_s liquid_apsk4; 1505 extern struct liquid_apsk_s liquid_apsk8; 1506 extern struct liquid_apsk_s liquid_apsk16; 1507 extern struct liquid_apsk_s liquid_apsk32; 1508 extern struct liquid_apsk_s liquid_apsk64; 1509 extern struct liquid_apsk_s liquid_apsk128; 1510 extern struct liquid_apsk_s liquid_apsk256; 1511 1512 1513 // 'square' 32-QAM (first quadrant) 1514 extern const float complex modem_arb_sqam32[8]; 1515 1516 // 'square' 128-QAM (first quadrant) 1517 extern const float complex modem_arb_sqam128[32]; 1518 1519 // V.29 star constellation 1520 extern const float complex modem_arb_V29[16]; 1521 1522 // Virginia Tech logo 1523 extern const float complex modem_arb_vt64[64]; 1524 1525 // optimal QAM constellations 1526 extern const float complex modem_arb16opt[16]; 1527 extern const float complex modem_arb32opt[32]; 1528 extern const float complex modem_arb64opt[64]; 1529 extern const float complex modem_arb128opt[128]; 1530 extern const float complex modem_arb256opt[256]; 1531 1532 1533 // 1534 // MODULE : multichannel 1535 // 1536 1537 // ofdm frame (common) 1538 1539 // generate short sequence symbols 1540 // _p : subcarrier allocation array 1541 // _M : total number of subcarriers 1542 // _S0 : output symbol (freq) 1543 // _s0 : output symbol (time) 1544 // _M_S0 : total number of enabled subcarriers in S0 1545 void ofdmframe_init_S0(unsigned char * _p, 1546 unsigned int _M, 1547 float complex * _S0, 1548 float complex * _s0, 1549 unsigned int * _M_S0); 1550 1551 // generate long sequence symbols 1552 // _p : subcarrier allocation array 1553 // _M : total number of subcarriers 1554 // _S1 : output symbol (freq) 1555 // _s1 : output symbol (time) 1556 // _M_S1 : total number of enabled subcarriers in S1 1557 void ofdmframe_init_S1(unsigned char * _p, 1558 unsigned int _M, 1559 float complex * _S1, 1560 float complex * _s1, 1561 unsigned int * _M_S1); 1562 1563 // generate symbol (add cyclic prefix/postfix, overlap) 1564 void ofdmframegen_gensymbol(ofdmframegen _q, 1565 float complex * _buffer); 1566 1567 void ofdmframesync_cpcorrelate(ofdmframesync _q); 1568 void ofdmframesync_findrxypeak(ofdmframesync _q); 1569 void ofdmframesync_rxpayload(ofdmframesync _q); 1570 1571 void ofdmframesync_execute_seekplcp(ofdmframesync _q); 1572 void ofdmframesync_execute_S0a(ofdmframesync _q); 1573 void ofdmframesync_execute_S0b(ofdmframesync _q); 1574 void ofdmframesync_execute_S1( ofdmframesync _q); 1575 void ofdmframesync_execute_rxsymbols(ofdmframesync _q); 1576 1577 void ofdmframesync_S0_metrics(ofdmframesync _q, 1578 float complex * _G, 1579 float complex * _s_hat); 1580 1581 // estimate short sequence gain 1582 // _q : ofdmframesync object 1583 // _x : input array (time) 1584 // _G : output gain (freq) 1585 void ofdmframesync_estimate_gain_S0(ofdmframesync _q, 1586 float complex * _x, 1587 float complex * _G); 1588 1589 // estimate long sequence gain 1590 // _q : ofdmframesync object 1591 // _x : input array (time) 1592 // _G : output gain (freq) 1593 void ofdmframesync_estimate_gain_S1(ofdmframesync _q, 1594 float complex * _x, 1595 float complex * _G); 1596 1597 // estimate complex equalizer gain from G0 and G1 1598 // _q : ofdmframesync object 1599 // _ntaps : number of time-domain taps for smoothing 1600 void ofdmframesync_estimate_eqgain(ofdmframesync _q, 1601 unsigned int _ntaps); 1602 1603 // estimate complex equalizer gain from G0 and G1 using polynomial fit 1604 // _q : ofdmframesync object 1605 // _order : polynomial order 1606 void ofdmframesync_estimate_eqgain_poly(ofdmframesync _q, 1607 unsigned int _order); 1608 1609 // recover symbol, correcting for gain, pilot phase, etc. 1610 void ofdmframesync_rxsymbol(ofdmframesync _q); 1611 1612 // 1613 // MODULE : nco (numerically-controlled oscillator) 1614 // 1615 1616 1617 // Numerically-controlled oscillator, floating point phase precision 1618 #define LIQUID_NCO_DEFINE_INTERNAL_API(NCO,T,TC) \ 1619 \ 1620 /* constrain phase/frequency to be in [-pi,pi) */ \ 1621 void NCO(_constrain_phase)(NCO() _q); \ 1622 void NCO(_constrain_frequency)(NCO() _q); \ 1623 \ 1624 /* compute trigonometric functions for nco/vco type */ \ 1625 void NCO(_compute_sincos_nco)(NCO() _q); \ 1626 void NCO(_compute_sincos_vco)(NCO() _q); \ 1627 \ 1628 /* reset internal phase-locked loop filter */ \ 1629 void NCO(_pll_reset)(NCO() _q); \ 1630 1631 // Define nco internal APIs 1632 LIQUID_NCO_DEFINE_INTERNAL_API(LIQUID_NCO_MANGLE_FLOAT, 1633 float, 1634 float complex) 1635 1636 // Numerically-controlled synthesizer (direct digital synthesis) 1637 #define LIQUID_SYNTH_DEFINE_INTERNAL_API(SYNTH,T,TC) \ 1638 \ 1639 /* constrain phase/frequency to be in [-pi,pi) */ \ 1640 void SYNTH(_constrain_phase)(SYNTH() _q); \ 1641 void SYNTH(_constrain_frequency)(SYNTH() _q); \ 1642 void SYNTH(_compute_synth)(SYNTH() _q); \ 1643 \ 1644 /* reset internal phase-locked loop filter */ \ 1645 void SYNTH(_pll_reset)(SYNTH() _q); \ 1646 1647 // Define nco internal APIs 1648 LIQUID_SYNTH_DEFINE_INTERNAL_API(SYNTH_MANGLE_FLOAT, 1649 float, 1650 liquid_float_complex) 1651 // 1652 // MODULE : optim (non-linear optimization) 1653 // 1654 1655 // optimization threshold switch 1656 // _u0 : first utility 1657 // _u1 : second utility 1658 // _minimize : minimize flag 1659 // 1660 // returns: 1661 // (_u0 > _u1) if (_minimize == 1) 1662 // (_u0 < _u1) otherwise 1663 int optim_threshold_switch(float _u0, 1664 float _u1, 1665 int _minimize); 1666 1667 // compute the gradient of a function at a particular point 1668 // _utility : user-defined function 1669 // _userdata : user-defined data object 1670 // _x : operating point, [size: _n x 1] 1671 // _n : dimensionality of search 1672 // _delta : step value for which to compute gradient 1673 // _gradient : resulting gradient 1674 void gradsearch_gradient(utility_function _utility, 1675 void * _userdata, 1676 float * _x, 1677 unsigned int _n, 1678 float _delta, 1679 float * _gradient); 1680 1681 // execute line search; loosely solve: 1682 // 1683 // min|max phi(alpha) := f(_x - alpha*_p) 1684 // 1685 // and return best guess at alpha that achieves this 1686 // 1687 // _utility : user-defined function 1688 // _userdata : user-defined data object 1689 // _direction : search direction (e.g. LIQUID_OPTIM_MINIMIZE) 1690 // _n : dimensionality of search 1691 // _x : operating point, [size: _n x 1] 1692 // _p : normalized gradient, [size: _n x 1] 1693 // _alpha : initial step size 1694 float gradsearch_linesearch(utility_function _utility, 1695 void * _userdata, 1696 int _direction, 1697 unsigned int _n, 1698 float * _x, 1699 float * _p, 1700 float _alpha); 1701 1702 // normalize vector, returning its l2-norm 1703 float gradsearch_norm(float * _v, 1704 unsigned int _n); 1705 1706 1707 // quasi-Newton search object 1708 struct qnsearch_s { 1709 float* v; // vector to optimize (externally allocated) 1710 unsigned int num_parameters; // number of parameters to optimize [n] 1711 1712 float gamma; // nominal stepsize 1713 float delta; // differential used to compute (estimate) derivative 1714 float dgamma; // decremental gamma parameter 1715 float gamma_hat; // step size (decreases each epoch) 1716 float* v_prime; // temporary vector array 1717 float* dv; // parameter step vector 1718 1719 float * B; // approximate Hessian matrix inverse [n x n] 1720 float * H; // Hessian matrix 1721 1722 float* p; // search direction 1723 float* gradient; // gradient approximation 1724 float* gradient0; // gradient approximation (previous step) 1725 1726 // External utility function. 1727 utility_function get_utility; 1728 float utility; // current utility 1729 void * userdata; // userdata pointer passed to utility callback 1730 int minimize; // minimize/maximimze utility (search direction) 1731 }; 1732 1733 // compute gradient(x_k) 1734 void qnsearch_compute_gradient(qnsearch _q); 1735 1736 // compute the norm of the gradient(x_k) 1737 void qnsearch_normalize_gradient(qnsearch _q); 1738 1739 // compute Hessian (estimate) 1740 void qnsearch_compute_Hessian(qnsearch _q); 1741 1742 // compute the updated inverse hessian matrix using the Broyden, Fletcher, 1743 // Goldfarb & Shanno method (BFGS) 1744 void qnsearch_update_hessian_bfgs(qnsearch _q); 1745 1746 1747 // Chromosome structure used in genetic algorithm searches 1748 struct chromosome_s { 1749 unsigned int num_traits; // number of represented traits 1750 unsigned int * bits_per_trait; // bits to represent each trait 1751 unsigned long * max_value; // maximum representable integer value 1752 unsigned long * traits; // chromosome data 1753 1754 unsigned int num_bits; // total number of bits 1755 }; 1756 1757 struct gasearch_s { 1758 chromosome * population; // population of chromosomes 1759 unsigned int population_size; // size of the population 1760 unsigned int selection_size; // number of 1761 float mutation_rate; // rate of mutation 1762 1763 unsigned int num_parameters; // number of parameters to optimize 1764 unsigned int bits_per_chromosome; // total number of bits in each chromosome 1765 1766 float *utility; // utility array 1767 unsigned int *rank; // rank indices of chromosomes (best to worst) 1768 1769 chromosome c; // copy of best chromosome, optimal solution 1770 float utility_opt; // optimum utility (fitness of best solution) 1771 1772 // External utility function. 1773 // 1774 // The success of a GA search algorithm is contingent upon the 1775 // design of a good utility function. It should meet the following 1776 // criteria: 1777 // - monotonically increasing (never flat) 1778 // - efficient to compute 1779 // - maps the [0,1] bounded output vector to desired range 1780 // - for multiple objectives, utility should be high \em only when 1781 // all objectives are met (multiplicative, not additive) 1782 gasearch_utility get_utility; // utility function pointer 1783 void * userdata; // object to optimize 1784 int minimize; // minimize/maximize utility (search direction) 1785 }; 1786 1787 // 1788 // gasearch internal methods 1789 // 1790 1791 // evaluate fitness of entire population 1792 void gasearch_evaluate(gasearch _q); 1793 1794 // crossover population 1795 void gasearch_crossover(gasearch _q); 1796 1797 // mutate population 1798 void gasearch_mutate(gasearch _q); 1799 1800 // rank population by fitness 1801 void gasearch_rank(gasearch _q); 1802 1803 // sort values by index 1804 // _v : input values [size: _len x 1] 1805 // _rank : output rank array (indices) [size: _len x 1] 1806 // _len : length of input array 1807 // _descending : descending/ascending 1808 void optim_sort(float *_v, 1809 unsigned int * _rank, 1810 unsigned int _len, 1811 int _descending); 1812 1813 1814 // 1815 // MODULE : random 1816 // 1817 1818 #define randf_inline() ((float) rand() / (float) RAND_MAX) 1819 1820 float complex icrandnf(); 1821 1822 // generate x ~ Gamma(delta,1) 1823 float randgammaf_delta(float _delta); 1824 1825 // data scrambler masks 1826 #define LIQUID_SCRAMBLE_MASK0 (0xb4) 1827 #define LIQUID_SCRAMBLE_MASK1 (0x6a) 1828 #define LIQUID_SCRAMBLE_MASK2 (0x8b) 1829 #define LIQUID_SCRAMBLE_MASK3 (0xc5) 1830 1831 // 1832 // MODULE : sequence 1833 // 1834 1835 // maximal-length sequence 1836 struct msequence_s { 1837 unsigned int m; // length generator polynomial, shift register 1838 unsigned int g; // generator polynomial 1839 unsigned int a; // initial shift register state, default: 1 1840 1841 unsigned int n; // length of sequence, n = (2^m)-1 1842 unsigned int v; // shift register 1843 unsigned int b; // return bit 1844 }; 1845 1846 // Default msequence generator objects 1847 extern struct msequence_s msequence_default[16]; 1848 1849 1850 // 1851 // MODULE : utility 1852 // 1853 1854 // number of ones in a byte 1855 // 0 0000 0000 : 0 1856 // 1 0000 0001 : 1 1857 // 2 0000 0010 : 1 1858 // 3 0000 0011 : 2 1859 // 4 0000 0100 : 1 1860 // ... 1861 // 126 0111 1110 : 6 1862 // 127 0111 1111 : 7 1863 // 128 1000 0000 : 1 1864 // 129 1000 0001 : 2 1865 // ... 1866 // 253 1111 1101 : 7 1867 // 254 1111 1110 : 7 1868 // 255 1111 1111 : 8 1869 extern const unsigned char liquid_c_ones[256]; 1870 1871 // Count the number of ones in an integer, inline insertion 1872 #define liquid_count_ones_uint16(x) ( \ 1873 liquid_c_ones[ (x) & 0xff ] + \ 1874 liquid_c_ones[ ((x)>>8) & 0xff ]) 1875 1876 #define liquid_count_ones_uint24(x) ( \ 1877 liquid_c_ones[ (x) & 0xff ] + \ 1878 liquid_c_ones[ ((x)>> 8) & 0xff ] + \ 1879 liquid_c_ones[ ((x)>>16) & 0xff ]) 1880 1881 #define liquid_count_ones_uint32(x) ( \ 1882 liquid_c_ones[ (x) & 0xff ] + \ 1883 liquid_c_ones[ ((x)>> 8) & 0xff ] + \ 1884 liquid_c_ones[ ((x)>>16) & 0xff ] + \ 1885 liquid_c_ones[ ((x)>>24) & 0xff ]) 1886 1887 1888 // number of ones in a byte, modulo 2 1889 // 0 0000 0000 : 0 1890 // 1 0000 0001 : 1 1891 // 2 0000 0010 : 1 1892 // 3 0000 0011 : 0 1893 // 4 0000 0100 : 1 1894 // ... 1895 // 126 0111 1110 : 0 1896 // 127 0111 1111 : 1 1897 // 128 1000 0000 : 1 1898 // 129 1000 0001 : 0 1899 // ... 1900 // 253 1111 1101 : 1 1901 // 254 1111 1110 : 1 1902 // 255 1111 1111 : 0 1903 extern const unsigned char liquid_c_ones_mod2[256]; 1904 1905 // Count the number of ones in an integer modulo 2, inline insertion 1906 #define liquid_count_ones_mod2_uint16(x) (( \ 1907 liquid_c_ones_mod2[ (x) & 0xff ] + \ 1908 liquid_c_ones_mod2[ ((x)>>8) & 0xff ]) % 2) 1909 1910 #define liquid_count_ones_mod2_uint32(x) (( \ 1911 liquid_c_ones_mod2[ (x) & 0xff ] + \ 1912 liquid_c_ones_mod2[ ((x)>> 8) & 0xff ] + \ 1913 liquid_c_ones_mod2[ ((x)>>16) & 0xff ] + \ 1914 liquid_c_ones_mod2[ ((x)>>24) & 0xff ]) % 2) 1915 1916 // compute binary dot-products (inline pre-processor macros) 1917 #define liquid_bdotprod_uint8(x,y) liquid_c_ones_mod2[(x)&(y)] 1918 #define liquid_bdotprod_uint16(x,y) liquid_count_ones_mod2_uint16((x)&(y)) 1919 #define liquid_bdotprod_uint32(x,y) liquid_count_ones_mod2_uint32((x)&(y)) 1920 1921 // number of leading zeros in byte 1922 extern unsigned int liquid_c_leading_zeros[256]; 1923 1924 // byte reversal and manipulation 1925 extern const unsigned char liquid_reverse_byte_gentab[256]; 1926 #endif // __LIQUID_INTERNAL_H__ 1927 1928