sys/tools/sound/feeder_rate_mkfilter.awk

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