1#!/usr/bin/awk -f 2# 3# SPDX-License-Identifier: BSD-2-Clause 4# 5# Copyright (c) 2007-2009 Ariff Abdullah <ariff@FreeBSD.org> 6# All rights reserved. 7# 8# Redistribution and use in source and binary forms, with or without 9# modification, are permitted provided that the following conditions 10# are met: 11# 1. Redistributions of source code must retain the above copyright 12# notice, this list of conditions and the following disclaimer. 13# 2. Redistributions in binary form must reproduce the above copyright 14# notice, this list of conditions and the following disclaimer in the 15# documentation and/or other materials provided with the distribution. 16# 17# THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS ``AS IS'' AND 18# ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE 19# IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE 20# ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE 21# FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL 22# DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS 23# OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) 24# HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT 25# LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY 26# OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF 27# SUCH DAMAGE. 28# 29# 30 31# 32# FIR filter design by windowing method. This might become one of the 33# funniest joke I've ever written due to too many tricks being applied to 34# ensure maximum precision (well, in fact this is already have the same 35# precision granularity compared to its C counterpart). Nevertheless, it's 36# working, precise, dynamically tunable based on "presets". 37# 38# XXX EXPECT TOTAL REWRITE! DON'T ARGUE! 39# 40# TODO: Using ultraspherical window might be a good idea. 41# 42# Based on: 43# 44# "Digital Audio Resampling" by Julius O. Smith III 45# 46# - http://ccrma.stanford.edu/~jos/resample/ 47# 48 49# 50# Some basic Math functions. 51# 52function abs(x) 53{ 54 return (((x < 0) ? -x : x) + 0); 55} 56 57function fabs(x) 58{ 59 return (((x < 0.0) ? -x : x) + 0.0); 60} 61 62function ceil(x, r) 63{ 64 r = int(x); 65 if (r < x) 66 r++; 67 return (r + 0); 68} 69 70function floor(x, r) 71{ 72 r = int(x); 73 if (r > x) 74 r--; 75 return (r + 0); 76} 77 78function pow(x, y) 79{ 80 return (exp(1.0 * y * log(1.0 * x))); 81} 82 83# 84# What the hell... 85# 86function shl(x, y) 87{ 88 while (y > 0) { 89 x *= 2; 90 y--; 91 } 92 return (x); 93} 94 95function shr(x, y) 96{ 97 while (y > 0 && x != 0) { 98 x = floor(x / 2); 99 y--; 100 } 101 return (x); 102} 103 104function fx_floor(v, o, r) 105{ 106 if (fabs(v) < fabs(smallest)) 107 smallest = v; 108 if (fabs(v) > fabs(largest)) 109 largest = v; 110 111 r = floor((v * o) + 0.5); 112 if (r < INT32_MIN || r > INT32_MAX) 113 printf("\n#error overflow v=%f, please reduce %d\n", v, o); 114 115 return (r); 116} 117 118# 119# Kaiser linear piecewise functions. 120# 121function kaiserAttn2Beta(attn, beta) 122{ 123 if (attn < 0.0) 124 return (Z_KAISER_BETA_DEFAULT); 125 126 if (attn > 50.0) 127 beta = 0.1102 * ((1.0 * attn) - 8.7); 128 else if (attn > 21.0) 129 beta = (0.5842 * pow((1.0 * attn) - 21.0, 0.4)) + \ 130 (0.07886 * ((1.0 * attn) - 21.0)); 131 else 132 beta = 0.0; 133 134 return (beta); 135} 136 137function kaiserBeta2Attn(beta, x, y, i, attn, xbeta) 138{ 139 if (beta < Z_WINDOW_KAISER) 140 return (Z_KAISER_ATTN_DEFAULT); 141 142 if (beta > kaiserAttn2Beta(50.0)) 143 attn = ((1.0 * beta) / 0.1102) + 8.7; 144 else { 145 x = 21.0; 146 y = 50.0; 147 attn = 0.5 * (x + y); 148 for (i = 0; i < 128; i++) { 149 xbeta = kaiserAttn2Beta(attn) 150 if (beta == xbeta || \ 151 (i > 63 && \ 152 fabs(beta - xbeta) < Z_KAISER_EPSILON)) 153 break; 154 if (beta > xbeta) 155 x = attn; 156 else 157 y = attn; 158 attn = 0.5 * (x + y); 159 } 160 } 161 162 return (attn); 163} 164 165function kaiserRolloff(len, attn) 166{ 167 return (1.0 - (((1.0 * attn) - 7.95) / (((1.0 * len) - 1.0) * 14.36))); 168} 169 170# 171# 0th order modified Bessel function of the first kind. 172# 173function I0(x, s, u, n, h, t) 174{ 175 s = n = u = 1.0; 176 h = x * 0.5; 177 178 do { 179 t = h / n; 180 n += 1.0; 181 t *= t; 182 u *= t; 183 s += u; 184 } while (u >= (I0_EPSILON * s)); 185 186 return (s); 187} 188 189function wname(beta) 190{ 191 if (beta >= Z_WINDOW_KAISER) 192 return ("Kaiser"); 193 else if (beta == Z_WINDOW_BLACKMAN_NUTTALL) 194 return ("Blackman - Nuttall"); 195 else if (beta == Z_WINDOW_NUTTALL) 196 return ("Nuttall"); 197 else if (beta == Z_WINDOW_BLACKMAN_HARRIS) 198 return ("Blackman - Harris"); 199 else if (beta == Z_WINDOW_BLACKMAN) 200 return ("Blackman"); 201 else if (beta == Z_WINDOW_HAMMING) 202 return ("Hamming"); 203 else if (beta == Z_WINDOW_HANN) 204 return ("Hann"); 205 else 206 return ("What The Hell !?!?"); 207} 208 209function rolloff_round(x) 210{ 211 if (x < 0.67) 212 x = 0.67; 213 else if (x > 1.0) 214 x = 1.0; 215 216 return (x); 217} 218 219function tap_round(x, y) 220{ 221 y = floor(x + 3); 222 y -= y % 4; 223 return (y); 224} 225 226function lpf(imp, n, rolloff, beta, num, i, j, x, nm, ibeta, w) 227{ 228 rolloff = rolloff_round(rolloff + (Z_NYQUIST_HOVER * (1.0 - rolloff))); 229 imp[0] = rolloff; 230 231 # 232 # Generate ideal sinc impulses, locate the last zero-crossing and pad 233 # the remaining with 0. 234 # 235 # Note that there are other (faster) ways of calculating this without 236 # the misery of traversing the entire sinc given the fact that the 237 # distance between each zero crossings is actually the bandwidth of 238 # the impulses, but it seems having 0.0001% chances of failure due to 239 # limited precision. 240 # 241 j = n; 242 for (i = 1; i < n; i++) { 243 x = (M_PI * i) / (1.0 * num); 244 imp[i] = sin(x * rolloff) / x; 245 if (i != 1 && (imp[i] * imp[i - 1]) <= 0.0) 246 j = i; 247 } 248 249 for (i = j; i < n; i++) 250 imp[i] = 0.0; 251 252 nm = 1.0 * (j - 1); 253 254 if (beta >= Z_WINDOW_KAISER) 255 ibeta = I0(beta); 256 257 for (i = 1; i < j; i++) { 258 if (beta >= Z_WINDOW_KAISER) { 259 # Kaiser window... 260 x = (1.0 * i) / nm; 261 w = I0(beta * sqrt(1.0 - (x * x))) / ibeta; 262 } else { 263 # Cosined windows... 264 x = (M_PI * i) / nm; 265 if (beta == Z_WINDOW_BLACKMAN_NUTTALL) { 266 # Blackman - Nuttall 267 w = 0.36335819 + (0.4891775 * cos(x)) + \ 268 (0.1365995 * cos(2 * x)) + \ 269 (0.0106411 * cos(3 * x)); 270 } else if (beta == Z_WINDOW_NUTTALL) { 271 # Nuttall 272 w = 0.355768 + (0.487396 * cos(x)) + \ 273 (0.144232 * cos(2 * x)) + \ 274 (0.012604 * cos(3 * x)); 275 } else if (beta == Z_WINDOW_BLACKMAN_HARRIS) { 276 # Blackman - Harris 277 w = 0.422323 + (0.49755 * cos(x)) + \ 278 (0.07922 * cos(2 * x)); 279 } else if (beta == Z_WINDOW_BLACKMAN) { 280 # Blackman 281 w = 0.42 + (0.50 * cos(x)) + \ 282 (0.08 * cos(2 * x)); 283 } else if (beta == Z_WINDOW_HAMMING) { 284 # Hamming 285 w = 0.54 + (0.46 * cos(x)); 286 } else if (beta == Z_WINDOW_HANN) { 287 # Hann 288 w = 0.50 + (0.50 * cos(x)); 289 } else { 290 # What The Hell !?!? 291 w = 0.0; 292 } 293 } 294 imp[i] *= w; 295 } 296 297 imp["impulse_length"] = j; 298 imp["rolloff"] = rolloff; 299} 300 301function mkfilter(imp, nmult, rolloff, beta, num, \ 302 nwing, mwing, nrolloff, i, dcgain, v, quality) 303{ 304 nwing = floor((nmult * num) / 2) + 1; 305 306 lpf(imp, nwing, rolloff, beta, num); 307 308 mwing = imp["impulse_length"]; 309 nrolloff = imp["rolloff"]; 310 quality = imp["quality"]; 311 312 dcgain = 0.0; 313 for (i = num; i < mwing; i += num) 314 dcgain += imp[i]; 315 dcgain *= 2.0; 316 dcgain += imp[0]; 317 318 for (i = 0; i < nwing; i++) 319 imp[i] /= dcgain; 320 321 if (quality > 2) 322 printf("\n"); 323 printf("/*\n"); 324 printf(" * quality = %d\n", quality); 325 printf(" * window = %s\n", wname(beta)); 326 if (beta >= Z_WINDOW_KAISER) { 327 printf(" * beta: %.2f\n", beta); 328 printf(" * stop: -%.2f dB\n", \ 329 kaiserBeta2Attn(beta)); 330 } 331 printf(" * length = %d\n", nmult); 332 printf(" * bandwidth = %.2f%%", rolloff * 100.0); 333 if (rolloff != nrolloff) { 334 printf(" + %.2f%% = %.2f%% (nyquist hover: %.2f%%)", \ 335 (nrolloff - rolloff) * 100.0, nrolloff * 100.0, \ 336 Z_NYQUIST_HOVER * 100.0); 337 } 338 printf("\n"); 339 printf(" * drift = %d\n", num); 340 printf(" * width = %d\n", mwing); 341 printf(" */\n"); 342 printf("static int32_t z_coeff_q%d[%d] = {", \ 343 quality, nwing + (Z_COEFF_OFFSET * 2)); 344 for (i = 0; i < (nwing + (Z_COEFF_OFFSET * 2)); i++) { 345 if ((i % 5) == 0) 346 printf("\n "); 347 if (i < Z_COEFF_OFFSET) 348 v = fx_floor(imp[Z_COEFF_OFFSET - i], Z_COEFF_ONE); 349 else if ((i - Z_COEFF_OFFSET) >= nwing) 350 v = fx_floor( \ 351 imp[nwing + nwing - i + Z_COEFF_OFFSET - 1],\ 352 Z_COEFF_ONE); 353 else 354 v = fx_floor(imp[i - Z_COEFF_OFFSET], Z_COEFF_ONE); 355 printf(" %s0x%08x,", (v < 0) ? "-" : " ", abs(v)); 356 } 357 printf("\n};\n\n"); 358 printf("/*\n"); 359 printf(" * interpolated q%d differences.\n", quality); 360 printf(" */\n"); 361 printf("static int32_t z_dcoeff_q%d[%d] = {", quality, nwing); 362 for (i = 1; i <= nwing; i++) { 363 if ((i % 5) == 1) 364 printf("\n "); 365 v = -imp[i - 1]; 366 if (i != nwing) 367 v += imp[i]; 368 v = fx_floor(v, Z_INTERP_COEFF_ONE); 369 if (abs(v) > abs(largest_interp)) 370 largest_interp = v; 371 printf(" %s0x%08x,", (v < 0) ? "-" : " ", abs(v)); 372 } 373 printf("\n};\n"); 374 375 return (nwing); 376} 377 378function filter_parse(s, a, i, attn, alen) 379{ 380 split(s, a, ":"); 381 alen = length(a); 382 383 if (alen > 0 && a[1] == "OVERSAMPLING_FACTOR") { 384 if (alen != 2) 385 return (-1); 386 init_drift(floor(a[2])); 387 return (-1); 388 } 389 390 if (alen > 0 && a[1] == "COEFFICIENT_BIT") { 391 if (alen != 2) 392 return (-1); 393 init_coeff_bit(floor(a[2])); 394 return (-1); 395 } 396 397 if (alen > 0 && a[1] == "ACCUMULATOR_BIT") { 398 if (alen != 2) 399 return (-1); 400 init_accum_bit(floor(a[2])); 401 return (-1); 402 } 403 404 if (alen > 0 && a[1] == "INTERPOLATOR") { 405 if (alen != 2) 406 return (-1); 407 init_coeff_interpolator(toupper(a[2])); 408 return (-1); 409 } 410 411 if (alen == 1 || alen == 2) { 412 if (a[1] == "NYQUIST_HOVER") { 413 i = 1.0 * a[2]; 414 Z_NYQUIST_HOVER = (i > 0.0 && i < 1.0) ? i : 0.0; 415 return (-1); 416 } 417 i = 1; 418 if (alen == 1) { 419 attn = Z_KAISER_ATTN_DEFAULT; 420 Popts["beta"] = Z_KAISER_BETA_DEFAULT; 421 } else if (Z_WINDOWS[a[1]] < Z_WINDOW_KAISER) { 422 Popts["beta"] = Z_WINDOWS[a[1]]; 423 i = tap_round(a[2]); 424 Popts["nmult"] = i; 425 if (i < 28) 426 i = 28; 427 i = 1.0 - (6.44 / i); 428 Popts["rolloff"] = rolloff_round(i); 429 return (0); 430 } else { 431 attn = 1.0 * a[i++]; 432 Popts["beta"] = kaiserAttn2Beta(attn); 433 } 434 i = tap_round(a[i]); 435 Popts["nmult"] = i; 436 if (i > 7 && i < 28) 437 i = 27; 438 i = kaiserRolloff(i, attn); 439 Popts["rolloff"] = rolloff_round(i); 440 441 return (0); 442 } 443 444 if (!(alen == 3 || alen == 4)) 445 return (-1); 446 447 i = 2; 448 449 if (a[1] == "kaiser") { 450 if (alen > 2) 451 Popts["beta"] = 1.0 * a[i++]; 452 else 453 Popts["beta"] = Z_KAISER_BETA_DEFAULT; 454 } else if (Z_WINDOWS[a[1]] < Z_WINDOW_KAISER) 455 Popts["beta"] = Z_WINDOWS[a[1]]; 456 else if (1.0 * a[1] < Z_WINDOW_KAISER) 457 return (-1); 458 else 459 Popts["beta"] = kaiserAttn2Beta(1.0 * a[1]); 460 Popts["nmult"] = tap_round(a[i++]); 461 if (a[1] == "kaiser" && alen == 3) 462 i = kaiserRolloff(Popts["nmult"], \ 463 kaiserBeta2Attn(Popts["beta"])); 464 else 465 i = 1.0 * a[i]; 466 Popts["rolloff"] = rolloff_round(i); 467 468 return (0); 469} 470 471function genscale(bit, s1, s2, scale) 472{ 473 if ((bit + Z_COEFF_SHIFT) > Z_ACCUMULATOR_BIT) 474 s1 = Z_COEFF_SHIFT - \ 475 (32 - (Z_ACCUMULATOR_BIT - Z_COEFF_SHIFT)); 476 else 477 s1 = Z_COEFF_SHIFT - (32 - bit); 478 479 s2 = Z_SHIFT + (32 - bit); 480 481 if (s1 == 0) 482 scale = "v"; 483 else if (s1 < 0) 484 scale = sprintf("(v) << %d", abs(s1)); 485 else 486 scale = sprintf("(v) >> %d", s1); 487 488 scale = sprintf("(%s) * Z_SCALE_CAST(s)", scale); 489 490 if (s2 != 0) 491 scale = sprintf("(%s) >> %d", scale, s2); 492 493 printf("#define Z_SCALE_%d(v, s)\t%s(%s)\n", \ 494 bit, (bit < 10) ? "\t" : "", scale); 495} 496 497function genlerp(bit, use64, lerp) 498{ 499 if ((bit + Z_LINEAR_SHIFT) <= 32) { 500 lerp = sprintf("(((y) - (x)) * (z)) >> %d", Z_LINEAR_SHIFT); 501 } else if (use64 != 0) { 502 if ((bit + Z_LINEAR_SHIFT) <= 64) { 503 lerp = sprintf( \ 504 "(((int64_t)(y) - (x)) * (z)) " \ 505 ">> %d", \ 506 Z_LINEAR_SHIFT); 507 } else { 508 lerp = sprintf( \ 509 "((int64_t)((y) >> %d) - ((x) >> %d)) * ", \ 510 "(z)" \ 511 bit + Z_LINEAR_SHIFT - 64, \ 512 bit + Z_LINEAR_SHIFT - 64); 513 if ((64 - bit) != 0) 514 lerp = sprintf("(%s) >> %d", lerp, 64 - bit); 515 } 516 } else { 517 lerp = sprintf( \ 518 "(((y) >> %d) - ((x) >> %d)) * (z)", \ 519 bit + Z_LINEAR_SHIFT - 32, \ 520 bit + Z_LINEAR_SHIFT - 32); 521 if ((32 - bit) != 0) 522 lerp = sprintf("(%s) >> %d", lerp, 32 - bit); 523 } 524 525 printf("#define Z_LINEAR_INTERPOLATE_%d(z, x, y)" \ 526 "\t\t\t\t%s\\\n\t((x) + (%s))\n", \ 527 bit, (bit < 10) ? "\t" : "", lerp); 528} 529 530function init_drift(drift, xdrift) 531{ 532 xdrift = floor(drift); 533 534 if (Z_DRIFT_SHIFT != -1) { 535 if (xdrift != Z_DRIFT_SHIFT) 536 printf("#error Z_DRIFT_SHIFT reinitialize!\n"); 537 return; 538 } 539 540 # 541 # Initialize filter oversampling factor, or in other word 542 # Z_DRIFT_SHIFT. 543 # 544 if (xdrift < 0) 545 xdrift = 1; 546 else if (xdrift > 31) 547 xdrift = 31; 548 549 Z_DRIFT_SHIFT = xdrift; 550 Z_DRIFT_ONE = shl(1, Z_DRIFT_SHIFT); 551 552 Z_SHIFT = Z_FULL_SHIFT - Z_DRIFT_SHIFT; 553 Z_ONE = shl(1, Z_SHIFT); 554 Z_MASK = Z_ONE - 1; 555} 556 557function init_coeff_bit(cbit, xcbit) 558{ 559 xcbit = floor(cbit); 560 561 if (Z_COEFF_SHIFT != 0) { 562 if (xcbit != Z_COEFF_SHIFT) 563 printf("#error Z_COEFF_SHIFT reinitialize!\n"); 564 return; 565 } 566 567 # 568 # Initialize dynamic range of coefficients. 569 # 570 if (xcbit < 1) 571 xcbit = 1; 572 else if (xcbit > 30) 573 xcbit = 30; 574 575 Z_COEFF_SHIFT = xcbit; 576 Z_COEFF_ONE = shl(1, Z_COEFF_SHIFT); 577} 578 579function init_accum_bit(accbit, xaccbit) 580{ 581 xaccbit = floor(accbit); 582 583 if (Z_ACCUMULATOR_BIT != 0) { 584 if (xaccbit != Z_ACCUMULATOR_BIT) 585 printf("#error Z_ACCUMULATOR_BIT reinitialize!\n"); 586 return; 587 } 588 589 # 590 # Initialize dynamic range of accumulator. 591 # 592 if (xaccbit > 64) 593 xaccbit = 64; 594 else if (xaccbit < 32) 595 xaccbit = 32; 596 597 Z_ACCUMULATOR_BIT = xaccbit; 598} 599 600function init_coeff_interpolator(interp) 601{ 602 # 603 # Validate interpolator type. 604 # 605 if (interp == "ZOH" || interp == "LINEAR" || \ 606 interp == "QUADRATIC" || interp == "HERMITE" || \ 607 interp == "BSPLINE" || interp == "OPT32X" || \ 608 interp == "OPT16X" || interp == "OPT8X" || \ 609 interp == "OPT4X" || interp == "OPT2X") 610 Z_COEFF_INTERPOLATOR = interp; 611} 612 613BEGIN { 614 I0_EPSILON = 1e-21; 615 M_PI = atan2(0.0, -1.0); 616 617 INT32_MAX = 1 + ((shl(1, 30) - 1) * 2); 618 INT32_MIN = -1 - INT32_MAX; 619 620 Z_COEFF_OFFSET = 5; 621 622 Z_ACCUMULATOR_BIT_DEFAULT = 58; 623 Z_ACCUMULATOR_BIT = 0; 624 625 Z_FULL_SHIFT = 30; 626 Z_FULL_ONE = shl(1, Z_FULL_SHIFT); 627 628 Z_COEFF_SHIFT_DEFAULT = 30; 629 Z_COEFF_SHIFT = 0; 630 Z_COEFF_ONE = 0; 631 632 Z_COEFF_INTERPOLATOR = 0; 633 634 Z_INTERP_COEFF_SHIFT = 24; 635 Z_INTERP_COEFF_ONE = shl(1, Z_INTERP_COEFF_SHIFT); 636 637 Z_LINEAR_FULL_SHIFT = Z_FULL_SHIFT; 638 Z_LINEAR_FULL_ONE = shl(1, Z_LINEAR_FULL_SHIFT); 639 Z_LINEAR_SHIFT = 8; 640 Z_LINEAR_UNSHIFT = Z_LINEAR_FULL_SHIFT - Z_LINEAR_SHIFT; 641 Z_LINEAR_ONE = shl(1, Z_LINEAR_SHIFT) 642 643 Z_DRIFT_SHIFT_DEFAULT = 5; 644 Z_DRIFT_SHIFT = -1; 645 # meehhhh... let it overflow... 646 #Z_SCALE_SHIFT = 31; 647 #Z_SCALE_ONE = shl(1, Z_SCALE_SHIFT); 648 649 Z_WINDOW_KAISER = 0.0; 650 Z_WINDOW_BLACKMAN_NUTTALL = -1.0; 651 Z_WINDOW_NUTTALL = -2.0; 652 Z_WINDOW_BLACKMAN_HARRIS = -3.0; 653 Z_WINDOW_BLACKMAN = -4.0; 654 Z_WINDOW_HAMMING = -5.0; 655 Z_WINDOW_HANN = -6.0; 656 657 Z_WINDOWS["blackman_nuttall"] = Z_WINDOW_BLACKMAN_NUTTALL; 658 Z_WINDOWS["nuttall"] = Z_WINDOW_NUTTALL; 659 Z_WINDOWS["blackman_harris"] = Z_WINDOW_BLACKMAN_HARRIS; 660 Z_WINDOWS["blackman"] = Z_WINDOW_BLACKMAN; 661 Z_WINDOWS["hamming"] = Z_WINDOW_HAMMING; 662 Z_WINDOWS["hann"] = Z_WINDOW_HANN; 663 664 Z_KAISER_2_BLACKMAN_BETA = 8.568611; 665 Z_KAISER_2_BLACKMAN_NUTTALL_BETA = 11.98; 666 667 Z_KAISER_ATTN_DEFAULT = 100; 668 Z_KAISER_BETA_DEFAULT = kaiserAttn2Beta(Z_KAISER_ATTN_DEFAULT); 669 670 Z_KAISER_EPSILON = 1e-21; 671 672 # 673 # This is practically a joke. 674 # 675 Z_NYQUIST_HOVER = 0.0; 676 677 smallest = 10.000000; 678 largest = 0.000010; 679 largest_interp = 0; 680 681 if (ARGC < 2) { 682 ARGC = 1; 683 ARGV[ARGC++] = "100:8:0.85"; 684 ARGV[ARGC++] = "100:36:0.92"; 685 ARGV[ARGC++] = "100:164:0.97"; 686 #ARGV[ARGC++] = "100:8"; 687 #ARGV[ARGC++] = "100:16"; 688 #ARGV[ARGC++] = "100:32:0.7929"; 689 #ARGV[ARGC++] = "100:64:0.8990"; 690 #ARGV[ARGC++] = "100:128:0.9499"; 691 } 692 693 printf("#ifndef _FEEDER_RATE_GEN_H_\n"); 694 printf("#define _FEEDER_RATE_GEN_H_\n\n"); 695 printf("/*\n"); 696 printf(" * Generated using feeder_rate_mkfilter.awk, heaven, wind and awesome.\n"); 697 printf(" *\n"); 698 printf(" * DO NOT EDIT!\n"); 699 printf(" */\n\n"); 700 printf("#define FEEDER_RATE_PRESETS\t\""); 701 for (i = 1; i < ARGC; i++) 702 printf("%s%s", (i == 1) ? "" : " ", ARGV[i]); 703 printf("\"\n\n"); 704 imp["quality"] = 2; 705 for (i = 1; i < ARGC; i++) { 706 if (filter_parse(ARGV[i]) == 0) { 707 beta = Popts["beta"]; 708 nmult = Popts["nmult"]; 709 rolloff = Popts["rolloff"]; 710 if (Z_DRIFT_SHIFT == -1) 711 init_drift(Z_DRIFT_SHIFT_DEFAULT); 712 if (Z_COEFF_SHIFT == 0) 713 init_coeff_bit(Z_COEFF_SHIFT_DEFAULT); 714 if (Z_ACCUMULATOR_BIT == 0) 715 init_accum_bit(Z_ACCUMULATOR_BIT_DEFAULT); 716 ztab[imp["quality"] - 2] = \ 717 mkfilter(imp, nmult, rolloff, beta, Z_DRIFT_ONE); 718 imp["quality"]++; 719 } 720 } 721 722 printf("\n"); 723 # 724 # XXX 725 # 726 #if (length(ztab) > 0) { 727 # j = 0; 728 # for (i = 0; i < length(ztab); i++) { 729 # if (ztab[i] > j) 730 # j = ztab[i]; 731 # } 732 # printf("static int32_t z_coeff_zero[%d] = {", j); 733 # for (i = 0; i < j; i++) { 734 # if ((i % 19) == 0) 735 # printf("\n"); 736 # printf(" 0,"); 737 # } 738 # printf("\n};\n\n"); 739 #} 740 # 741 # XXX 742 # 743 printf("static const struct {\n"); 744 printf("\tint32_t len;\n"); 745 printf("\tint32_t *coeff;\n"); 746 printf("\tint32_t *dcoeff;\n"); 747 printf("} z_coeff_tab[] = {\n"); 748 if (length(ztab) > 0) { 749 j = 0; 750 for (i = 0; i < length(ztab); i++) { 751 if (ztab[i] > j) 752 j = ztab[i]; 753 } 754 j = length(sprintf("%d", j)); 755 lfmt = sprintf("%%%dd", j); 756 j = length(sprintf("z_coeff_q%d", length(ztab) + 1)); 757 zcfmt = sprintf("%%-%ds", j); 758 zdcfmt = sprintf("%%-%ds", j + 1); 759 760 for (i = 0; i < length(ztab); i++) { 761 l = sprintf(lfmt, ztab[i]); 762 zc = sprintf("z_coeff_q%d", i + 2); 763 zc = sprintf(zcfmt, zc); 764 zdc = sprintf("z_dcoeff_q%d", i + 2); 765 zdc = sprintf(zdcfmt, zdc); 766 printf("\t{ %s, %s, %s },\n", l, zc, zdc); 767 } 768 } else 769 printf("\t{ 0, NULL, NULL }\n"); 770 printf("};\n\n"); 771 772 #Z_UNSHIFT = 0; 773 #v = shr(Z_ONE - 1, Z_UNSHIFT) * abs(largest_interp); 774 #while (v < 0 || v > INT32_MAX) { 775 # Z_UNSHIFT += 1; 776 # v = shr(Z_ONE - 1, Z_UNSHIFT) * abs(largest_interp); 777 #} 778 v = ((Z_ONE - 1) * abs(largest_interp)) / INT32_MAX; 779 Z_UNSHIFT = ceil(log(v) / log(2.0)); 780 Z_INTERP_SHIFT = Z_SHIFT - Z_UNSHIFT + Z_INTERP_COEFF_SHIFT; 781 782 Z_INTERP_UNSHIFT = (Z_SHIFT - Z_UNSHIFT) + Z_INTERP_COEFF_SHIFT \ 783 - Z_COEFF_SHIFT; 784 785 printf("#define Z_COEFF_TAB_SIZE\t\t\t\t\t\t\\\n"); 786 printf("\t((int32_t)(sizeof(z_coeff_tab) /"); 787 printf(" sizeof(z_coeff_tab[0])))\n\n"); 788 printf("#define Z_COEFF_OFFSET\t\t%d\n\n", Z_COEFF_OFFSET); 789 printf("#define Z_RSHIFT(x, y)\t\t(((x) + " \ 790 "(1 << ((y) - 1))) >> (y))\n"); 791 printf("#define Z_RSHIFT_L(x, y)\t(((x) + " \ 792 "(1LL << ((y) - 1))) >> (y))\n\n"); 793 printf("#define Z_FULL_SHIFT\t\t%d\n", Z_FULL_SHIFT); 794 printf("#define Z_FULL_ONE\t\t0x%08x%s\n", Z_FULL_ONE, \ 795 (Z_FULL_ONE > INT32_MAX) ? "U" : ""); 796 printf("\n"); 797 printf("#define Z_DRIFT_SHIFT\t\t%d\n", Z_DRIFT_SHIFT); 798 #printf("#define Z_DRIFT_ONE\t\t0x%08x\n", Z_DRIFT_ONE); 799 printf("\n"); 800 printf("#define Z_SHIFT\t\t\t%d\n", Z_SHIFT); 801 printf("#define Z_ONE\t\t\t0x%08x\n", Z_ONE); 802 printf("#define Z_MASK\t\t\t0x%08x\n", Z_MASK); 803 printf("\n"); 804 printf("#define Z_COEFF_SHIFT\t\t%d\n", Z_COEFF_SHIFT); 805 zinterphp = "(z) * (d)"; 806 zinterpunshift = Z_SHIFT + Z_INTERP_COEFF_SHIFT - Z_COEFF_SHIFT; 807 if (zinterpunshift > 0) { 808 v = (Z_ONE - 1) * abs(largest_interp); 809 if (v < INT32_MIN || v > INT32_MAX) 810 zinterphp = sprintf("(int64_t)%s", zinterphp); 811 zinterphp = sprintf("(%s) >> %d", zinterphp, zinterpunshift); 812 } else if (zinterpunshift < 0) 813 zinterphp = sprintf("(%s) << %d", zinterphp, \ 814 abs(zinterpunshift)); 815 if (Z_UNSHIFT == 0) 816 zinterp = "z"; 817 else 818 zinterp = sprintf("(z) >> %d", Z_UNSHIFT); 819 zinterp = sprintf("(%s) * (d)", zinterp); 820 if (Z_INTERP_UNSHIFT < 0) 821 zinterp = sprintf("(%s) << %d", zinterp, \ 822 abs(Z_INTERP_UNSHIFT)); 823 else if (Z_INTERP_UNSHIFT > 0) 824 zinterp = sprintf("(%s) >> %d", zinterp, Z_INTERP_UNSHIFT); 825 if (zinterphp != zinterp) { 826 printf("\n#ifdef SND_FEEDER_RATE_HP\n"); 827 printf("#define Z_COEFF_INTERPOLATE(z, c, d)" \ 828 "\t\t\t\t\t\\\n\t((c) + (%s))\n", zinterphp); 829 printf("#else\n"); 830 printf("#define Z_COEFF_INTERPOLATE(z, c, d)" \ 831 "\t\t\t\t\t\\\n\t((c) + (%s))\n", zinterp); 832 printf("#endif\n"); 833 } else 834 printf("#define Z_COEFF_INTERPOLATE(z, c, d)" \ 835 "\t\t\t\t\t\\\n\t((c) + (%s))\n", zinterp); 836 #printf("\n"); 837 #printf("#define Z_SCALE_SHIFT\t\t%d\n", Z_SCALE_SHIFT); 838 #printf("#define Z_SCALE_ONE\t\t0x%08x%s\n", Z_SCALE_ONE, \ 839 # (Z_SCALE_ONE > INT32_MAX) ? "U" : ""); 840 printf("\n"); 841 printf("#define Z_SCALE_CAST(s)\t\t((uint32_t)(s))\n"); 842 genscale(8); 843 genscale(16); 844 genscale(24); 845 genscale(32); 846 printf("\n"); 847 printf("#define Z_ACCUMULATOR_BIT\t%d\n\n", Z_ACCUMULATOR_BIT) 848 for (i = 8; i <= 32; i += 8) { 849 gbit = ((i + Z_COEFF_SHIFT) > Z_ACCUMULATOR_BIT) ? \ 850 (i - (Z_ACCUMULATOR_BIT - Z_COEFF_SHIFT)) : 0; 851 printf("#define Z_GUARD_BIT_%d\t\t%d\n", i, gbit); 852 if (gbit == 0) 853 printf("#define Z_NORM_%d(v)\t\t(v)\n\n", i); 854 else 855 printf("#define Z_NORM_%d(v)\t\t" \ 856 "((v) >> Z_GUARD_BIT_%d)\n\n", i, i); 857 } 858 printf("\n"); 859 printf("#define Z_LINEAR_FULL_ONE\t0x%08xU\n", Z_LINEAR_FULL_ONE); 860 printf("#define Z_LINEAR_SHIFT\t\t%d\n", Z_LINEAR_SHIFT); 861 printf("#define Z_LINEAR_UNSHIFT\t%d\n", Z_LINEAR_UNSHIFT); 862 printf("#define Z_LINEAR_ONE\t\t0x%08x\n", Z_LINEAR_ONE); 863 printf("\n"); 864 printf("#ifdef SND_PCM_64\n"); 865 genlerp(8, 1); 866 genlerp(16, 1); 867 genlerp(24, 1); 868 genlerp(32, 1); 869 printf("#else\t/* !SND_PCM_64 */\n"); 870 genlerp(8, 0); 871 genlerp(16, 0); 872 genlerp(24, 0); 873 genlerp(32, 0); 874 printf("#endif\t/* SND_PCM_64 */\n"); 875 printf("\n"); 876 printf("#define Z_QUALITY_ZOH\t\t0\n"); 877 printf("#define Z_QUALITY_LINEAR\t1\n"); 878 printf("#define Z_QUALITY_SINC\t\t%d\n", \ 879 floor((length(ztab) - 1) / 2) + 2); 880 printf("\n"); 881 printf("#define Z_QUALITY_MIN\t\t0\n"); 882 printf("#define Z_QUALITY_MAX\t\t%d\n", length(ztab) + 1); 883 if (Z_COEFF_INTERPOLATOR != 0) 884 printf("\n#define Z_COEFF_INTERP_%s\t1\n", \ 885 Z_COEFF_INTERPOLATOR); 886 printf("\n/*\n * smallest: %.32f\n * largest: %.32f\n *\n", \ 887 smallest, largest); 888 printf(" * z_unshift=%d, z_interp_shift=%d\n *\n", \ 889 Z_UNSHIFT, Z_INTERP_SHIFT); 890 v = shr(Z_ONE - 1, Z_UNSHIFT) * abs(largest_interp); 891 printf(" * largest interpolation multiplication: %d\n */\n", v); 892 if (v < INT32_MIN || v > INT32_MAX) { 893 printf("\n#ifndef SND_FEEDER_RATE_HP\n"); 894 printf("#error interpolation overflow, please reduce" \ 895 " Z_INTERP_SHIFT\n"); 896 printf("#endif\n"); 897 } 898 899 printf("\n#endif /* !_FEEDER_RATE_GEN_H_ */\n"); 900} 901