1193640Sariff#!/usr/bin/awk -f 2193640Sariff# 3193640Sariff# Copyright (c) 2007-2009 Ariff Abdullah <ariff@FreeBSD.org> 4193640Sariff# All rights reserved. 5193640Sariff# 6193640Sariff# Redistribution and use in source and binary forms, with or without 7193640Sariff# modification, are permitted provided that the following conditions 8193640Sariff# are met: 9193640Sariff# 1. Redistributions of source code must retain the above copyright 10193640Sariff# notice, this list of conditions and the following disclaimer. 11193640Sariff# 2. Redistributions in binary form must reproduce the above copyright 12193640Sariff# notice, this list of conditions and the following disclaimer in the 13193640Sariff# documentation and/or other materials provided with the distribution. 14193640Sariff# 15193640Sariff# THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS ``AS IS'' AND 16193640Sariff# ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE 17193640Sariff# IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE 18193640Sariff# ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE 19193640Sariff# FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL 20193640Sariff# DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS 21193640Sariff# OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) 22193640Sariff# HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT 23193640Sariff# LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY 24193640Sariff# OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF 25193640Sariff# SUCH DAMAGE. 26193640Sariff# 27193640Sariff# $FreeBSD$ 28193640Sariff# 29193640Sariff 30193640Sariff# 31193640Sariff# FIR filter design by windowing method. This might become one of the 32193640Sariff# funniest joke I've ever written due to too many tricks being applied to 33193640Sariff# ensure maximum precision (well, in fact this is already have the same 34193640Sariff# precision granularity compared to its C counterpart). Nevertheless, it's 35193640Sariff# working, precise, dynamically tunable based on "presets". 36193640Sariff# 37193640Sariff# XXX EXPECT TOTAL REWRITE! DON'T ARGUE! 38193640Sariff# 39193640Sariff# TODO: Using ultraspherical window might be a good idea. 40193640Sariff# 41193640Sariff# Based on: 42193640Sariff# 43193640Sariff# "Digital Audio Resampling" by Julius O. Smith III 44193640Sariff# 45193640Sariff# - http://ccrma.stanford.edu/~jos/resample/ 46193640Sariff# 47193640Sariff 48193640Sariff# 49193640Sariff# Some basic Math functions. 50193640Sariff# 51193640Sarifffunction abs(x) 52193640Sariff{ 53193640Sariff return (((x < 0) ? -x : x) + 0); 54193640Sariff} 55193640Sariff 56193640Sarifffunction fabs(x) 57193640Sariff{ 58193640Sariff return (((x < 0.0) ? -x : x) + 0.0); 59193640Sariff} 60193640Sariff 61193640Sarifffunction ceil(x, r) 62193640Sariff{ 63193640Sariff r = int(x); 64193640Sariff if (r < x) 65193640Sariff r++; 66193640Sariff return (r + 0); 67193640Sariff} 68193640Sariff 69193640Sarifffunction floor(x, r) 70193640Sariff{ 71193640Sariff r = int(x); 72193640Sariff if (r > x) 73193640Sariff r--; 74193640Sariff return (r + 0); 75193640Sariff} 76193640Sariff 77193640Sarifffunction pow(x, y) 78193640Sariff{ 79193640Sariff return (exp(1.0 * y * log(1.0 * x))); 80193640Sariff} 81193640Sariff 82193640Sariff# 83193640Sariff# What the hell... 84193640Sariff# 85193640Sarifffunction shl(x, y) 86193640Sariff{ 87193640Sariff while (y > 0) { 88193640Sariff x *= 2; 89193640Sariff y--; 90193640Sariff } 91193640Sariff return (x); 92193640Sariff} 93193640Sariff 94193640Sarifffunction shr(x, y) 95193640Sariff{ 96193640Sariff while (y > 0 && x != 0) { 97193640Sariff x = floor(x / 2); 98193640Sariff y--; 99193640Sariff } 100193640Sariff return (x); 101193640Sariff} 102193640Sariff 103193640Sarifffunction fx_floor(v, o, r) 104193640Sariff{ 105193640Sariff if (fabs(v) < fabs(smallest)) 106193640Sariff smallest = v; 107193640Sariff if (fabs(v) > fabs(largest)) 108193640Sariff largest = v; 109193640Sariff 110193640Sariff r = floor((v * o) + 0.5); 111193640Sariff if (r < INT32_MIN || r > INT32_MAX) 112193640Sariff printf("\n#error overflow v=%f, please reduce %d\n", v, o); 113193640Sariff 114193640Sariff return (r); 115193640Sariff} 116193640Sariff 117193640Sariff# 118193640Sariff# Kaiser linear piecewise functions. 119193640Sariff# 120193640Sarifffunction kaiserAttn2Beta(attn, beta) 121193640Sariff{ 122193640Sariff if (attn < 0.0) 123193640Sariff return (Z_KAISER_BETA_DEFAULT); 124193640Sariff 125193640Sariff if (attn > 50.0) 126193640Sariff beta = 0.1102 * ((1.0 * attn) - 8.7); 127193640Sariff else if (attn > 21.0) 128193640Sariff beta = (0.5842 * pow((1.0 * attn) - 21.0, 0.4)) + \ 129193640Sariff (0.07886 * ((1.0 * attn) - 21.0)); 130193640Sariff else 131193640Sariff beta = 0.0; 132193640Sariff 133193640Sariff return (beta); 134193640Sariff} 135193640Sariff 136193640Sarifffunction kaiserBeta2Attn(beta, x, y, i, attn, xbeta) 137193640Sariff{ 138193640Sariff if (beta < Z_WINDOW_KAISER) 139193640Sariff return (Z_KAISER_ATTN_DEFAULT); 140193640Sariff 141193640Sariff if (beta > kaiserAttn2Beta(50.0)) 142193640Sariff attn = ((1.0 * beta) / 0.1102) + 8.7; 143193640Sariff else { 144193640Sariff x = 21.0; 145193640Sariff y = 50.0; 146193640Sariff attn = 0.5 * (x + y); 147193640Sariff for (i = 0; i < 128; i++) { 148193640Sariff xbeta = kaiserAttn2Beta(attn) 149193640Sariff if (beta == xbeta || \ 150193640Sariff (i > 63 && \ 151193640Sariff fabs(beta - xbeta) < Z_KAISER_EPSILON)) 152193640Sariff break; 153193640Sariff if (beta > xbeta) 154193640Sariff x = attn; 155193640Sariff else 156193640Sariff y = attn; 157193640Sariff attn = 0.5 * (x + y); 158193640Sariff } 159193640Sariff } 160193640Sariff 161193640Sariff return (attn); 162193640Sariff} 163193640Sariff 164193640Sarifffunction kaiserRolloff(len, attn) 165193640Sariff{ 166193640Sariff return (1.0 - (((1.0 * attn) - 7.95) / (((1.0 * len) - 1.0) * 14.36))); 167193640Sariff} 168193640Sariff 169193640Sariff# 170193640Sariff# 0th order modified Bessel function of the first kind. 171193640Sariff# 172193640Sarifffunction I0(x, s, u, n, h, t) 173193640Sariff{ 174193640Sariff s = n = u = 1.0; 175193640Sariff h = x * 0.5; 176193640Sariff 177193640Sariff do { 178193640Sariff t = h / n; 179193640Sariff n += 1.0; 180193640Sariff t *= t; 181193640Sariff u *= t; 182193640Sariff s += u; 183193640Sariff } while (u >= (I0_EPSILON * s)); 184193640Sariff 185193640Sariff return (s); 186193640Sariff} 187193640Sariff 188193640Sarifffunction wname(beta) 189193640Sariff{ 190193640Sariff if (beta >= Z_WINDOW_KAISER) 191193640Sariff return ("Kaiser"); 192193640Sariff else if (beta == Z_WINDOW_BLACKMAN_NUTTALL) 193193640Sariff return ("Blackman - Nuttall"); 194193640Sariff else if (beta == Z_WINDOW_NUTTALL) 195193640Sariff return ("Nuttall"); 196193640Sariff else if (beta == Z_WINDOW_BLACKMAN_HARRIS) 197193640Sariff return ("Blackman - Harris"); 198193640Sariff else if (beta == Z_WINDOW_BLACKMAN) 199193640Sariff return ("Blackman"); 200193640Sariff else if (beta == Z_WINDOW_HAMMING) 201193640Sariff return ("Hamming"); 202193640Sariff else if (beta == Z_WINDOW_HANN) 203193640Sariff return ("Hann"); 204193640Sariff else 205193640Sariff return ("What The Hell !?!?"); 206193640Sariff} 207193640Sariff 208193640Sarifffunction rolloff_round(x) 209193640Sariff{ 210193640Sariff if (x < 0.67) 211193640Sariff x = 0.67; 212193640Sariff else if (x > 1.0) 213193640Sariff x = 1.0; 214193640Sariff 215193640Sariff return (x); 216193640Sariff} 217193640Sariff 218193640Sarifffunction tap_round(x, y) 219193640Sariff{ 220193640Sariff y = floor(x + 3); 221193640Sariff y -= y % 4; 222193640Sariff return (y); 223193640Sariff} 224193640Sariff 225193640Sarifffunction lpf(imp, n, rolloff, beta, num, i, j, x, nm, ibeta, w) 226193640Sariff{ 227193640Sariff rolloff = rolloff_round(rolloff + (Z_NYQUIST_HOVER * (1.0 - rolloff))); 228193640Sariff imp[0] = rolloff; 229193640Sariff 230193640Sariff # 231193640Sariff # Generate ideal sinc impulses, locate the last zero-crossing and pad 232193640Sariff # the remaining with 0. 233193640Sariff # 234193640Sariff # Note that there are other (faster) ways of calculating this without 235193640Sariff # the misery of traversing the entire sinc given the fact that the 236193640Sariff # distance between each zero crossings is actually the bandwidth of 237193640Sariff # the impulses, but it seems having 0.0001% chances of failure due to 238193640Sariff # limited precision. 239193640Sariff # 240193640Sariff j = n; 241193640Sariff for (i = 1; i < n; i++) { 242193640Sariff x = (M_PI * i) / (1.0 * num); 243193640Sariff imp[i] = sin(x * rolloff) / x; 244193640Sariff if (i != 1 && (imp[i] * imp[i - 1]) <= 0.0) 245193640Sariff j = i; 246193640Sariff } 247193640Sariff 248193640Sariff for (i = j; i < n; i++) 249193640Sariff imp[i] = 0.0; 250193640Sariff 251193640Sariff nm = 1.0 * (j - 1); 252193640Sariff 253193640Sariff if (beta >= Z_WINDOW_KAISER) 254193640Sariff ibeta = I0(beta); 255193640Sariff 256193640Sariff for (i = 1; i < j; i++) { 257193640Sariff if (beta >= Z_WINDOW_KAISER) { 258193640Sariff # Kaiser window... 259193640Sariff x = (1.0 * i) / nm; 260193640Sariff w = I0(beta * sqrt(1.0 - (x * x))) / ibeta; 261193640Sariff } else { 262193640Sariff # Cosined windows... 263193640Sariff x = (M_PI * i) / nm; 264193640Sariff if (beta == Z_WINDOW_BLACKMAN_NUTTALL) { 265193640Sariff # Blackman - Nuttall 266193640Sariff w = 0.36335819 + (0.4891775 * cos(x)) + \ 267193640Sariff (0.1365995 * cos(2 * x)) + \ 268193640Sariff (0.0106411 * cos(3 * x)); 269193640Sariff } else if (beta == Z_WINDOW_NUTTALL) { 270193640Sariff # Nuttall 271193640Sariff w = 0.355768 + (0.487396 * cos(x)) + \ 272193640Sariff (0.144232 * cos(2 * x)) + \ 273193640Sariff (0.012604 * cos(3 * x)); 274193640Sariff } else if (beta == Z_WINDOW_BLACKMAN_HARRIS) { 275193640Sariff # Blackman - Harris 276193640Sariff w = 0.422323 + (0.49755 * cos(x)) + \ 277193640Sariff (0.07922 * cos(2 * x)); 278193640Sariff } else if (beta == Z_WINDOW_BLACKMAN) { 279193640Sariff # Blackman 280193640Sariff w = 0.42 + (0.50 * cos(x)) + \ 281193640Sariff (0.08 * cos(2 * x)); 282193640Sariff } else if (beta == Z_WINDOW_HAMMING) { 283193640Sariff # Hamming 284193640Sariff w = 0.54 + (0.46 * cos(x)); 285193640Sariff } else if (beta == Z_WINDOW_HANN) { 286193640Sariff # Hann 287193640Sariff w = 0.50 + (0.50 * cos(x)); 288193640Sariff } else { 289193640Sariff # What The Hell !?!? 290193640Sariff w = 0.0; 291193640Sariff } 292193640Sariff } 293193640Sariff imp[i] *= w; 294193640Sariff } 295193640Sariff 296193640Sariff imp["impulse_length"] = j; 297193640Sariff imp["rolloff"] = rolloff; 298193640Sariff} 299193640Sariff 300193640Sarifffunction mkfilter(imp, nmult, rolloff, beta, num, \ 301193640Sariff nwing, mwing, nrolloff, i, dcgain, v, quality) 302193640Sariff{ 303193640Sariff nwing = floor((nmult * num) / 2) + 1; 304193640Sariff 305193640Sariff lpf(imp, nwing, rolloff, beta, num); 306193640Sariff 307193640Sariff mwing = imp["impulse_length"]; 308193640Sariff nrolloff = imp["rolloff"]; 309193640Sariff quality = imp["quality"]; 310193640Sariff 311193640Sariff dcgain = 0.0; 312193640Sariff for (i = num; i < mwing; i += num) 313193640Sariff dcgain += imp[i]; 314193640Sariff dcgain *= 2.0; 315193640Sariff dcgain += imp[0]; 316193640Sariff 317193640Sariff for (i = 0; i < nwing; i++) 318193640Sariff imp[i] /= dcgain; 319193640Sariff 320193640Sariff if (quality > 2) 321193640Sariff printf("\n"); 322193640Sariff printf("/*\n"); 323193640Sariff printf(" * quality = %d\n", quality); 324193640Sariff printf(" * window = %s\n", wname(beta)); 325193640Sariff if (beta >= Z_WINDOW_KAISER) { 326193640Sariff printf(" * beta: %.2f\n", beta); 327193640Sariff printf(" * stop: -%.2f dB\n", \ 328193640Sariff kaiserBeta2Attn(beta)); 329193640Sariff } 330193640Sariff printf(" * length = %d\n", nmult); 331193640Sariff printf(" * bandwidth = %.2f%%", rolloff * 100.0); 332193640Sariff if (rolloff != nrolloff) { 333193640Sariff printf(" + %.2f%% = %.2f%% (nyquist hover: %.2f%%)", \ 334193640Sariff (nrolloff - rolloff) * 100.0, nrolloff * 100.0, \ 335193640Sariff Z_NYQUIST_HOVER * 100.0); 336193640Sariff } 337193640Sariff printf("\n"); 338193640Sariff printf(" * drift = %d\n", num); 339193640Sariff printf(" * width = %d\n", mwing); 340193640Sariff printf(" */\n"); 341193640Sariff printf("static int32_t z_coeff_q%d[%d] = {", \ 342193640Sariff quality, nwing + (Z_COEFF_OFFSET * 2)); 343193640Sariff for (i = 0; i < (nwing + (Z_COEFF_OFFSET * 2)); i++) { 344193640Sariff if ((i % 5) == 0) 345193640Sariff printf("\n "); 346193640Sariff if (i < Z_COEFF_OFFSET) 347193640Sariff v = fx_floor(imp[Z_COEFF_OFFSET - i], Z_COEFF_ONE); 348193640Sariff else if ((i - Z_COEFF_OFFSET) >= nwing) 349193640Sariff v = fx_floor( \ 350193640Sariff imp[nwing + nwing - i + Z_COEFF_OFFSET - 1],\ 351193640Sariff Z_COEFF_ONE); 352193640Sariff else 353193640Sariff v = fx_floor(imp[i - Z_COEFF_OFFSET], Z_COEFF_ONE); 354193640Sariff printf(" %s0x%08x,", (v < 0) ? "-" : " ", abs(v)); 355193640Sariff } 356193640Sariff printf("\n};\n\n"); 357193640Sariff printf("/*\n"); 358193640Sariff printf(" * interpolated q%d differences.\n", quality); 359193640Sariff printf(" */\n"); 360193640Sariff printf("static int32_t z_dcoeff_q%d[%d] = {", quality, nwing); 361193640Sariff for (i = 1; i <= nwing; i++) { 362193640Sariff if ((i % 5) == 1) 363193640Sariff printf("\n "); 364193640Sariff v = -imp[i - 1]; 365193640Sariff if (i != nwing) 366193640Sariff v += imp[i]; 367193640Sariff v = fx_floor(v, Z_INTERP_COEFF_ONE); 368193640Sariff if (abs(v) > abs(largest_interp)) 369193640Sariff largest_interp = v; 370193640Sariff printf(" %s0x%08x,", (v < 0) ? "-" : " ", abs(v)); 371193640Sariff } 372193640Sariff printf("\n};\n"); 373193640Sariff 374193640Sariff return (nwing); 375193640Sariff} 376193640Sariff 377193640Sarifffunction filter_parse(s, a, i, attn, alen) 378193640Sariff{ 379193640Sariff split(s, a, ":"); 380193640Sariff alen = length(a); 381193640Sariff 382194233Sariff if (alen > 0 && a[1] == "OVERSAMPLING_FACTOR") { 383194233Sariff if (alen != 2) 384194233Sariff return (-1); 385194233Sariff init_drift(floor(a[2])); 386194233Sariff return (-1); 387194233Sariff } 388194233Sariff 389195378Sariff if (alen > 0 && a[1] == "COEFFICIENT_BIT") { 390195378Sariff if (alen != 2) 391195378Sariff return (-1); 392195378Sariff init_coeff_bit(floor(a[2])); 393195378Sariff return (-1); 394195378Sariff } 395195378Sariff 396195378Sariff if (alen > 0 && a[1] == "ACCUMULATOR_BIT") { 397195378Sariff if (alen != 2) 398195378Sariff return (-1); 399195378Sariff init_accum_bit(floor(a[2])); 400195378Sariff return (-1); 401195378Sariff } 402195378Sariff 403195378Sariff if (alen > 0 && a[1] == "INTERPOLATOR") { 404195378Sariff if (alen != 2) 405195378Sariff return (-1); 406195378Sariff init_coeff_interpolator(toupper(a[2])); 407195378Sariff return (-1); 408195378Sariff } 409195378Sariff 410193640Sariff if (alen == 1 || alen == 2) { 411194233Sariff if (a[1] == "NYQUIST_HOVER") { 412193640Sariff i = 1.0 * a[2]; 413193640Sariff Z_NYQUIST_HOVER = (i > 0.0 && i < 1.0) ? i : 0.0; 414193640Sariff return (-1); 415193640Sariff } 416193640Sariff i = 1; 417193640Sariff if (alen == 1) { 418193640Sariff attn = Z_KAISER_ATTN_DEFAULT; 419193640Sariff Popts["beta"] = Z_KAISER_BETA_DEFAULT; 420193640Sariff } else if (Z_WINDOWS[a[1]] < Z_WINDOW_KAISER) { 421193640Sariff Popts["beta"] = Z_WINDOWS[a[1]]; 422193640Sariff i = tap_round(a[2]); 423193640Sariff Popts["nmult"] = i; 424193640Sariff if (i < 28) 425193640Sariff i = 28; 426193640Sariff i = 1.0 - (6.44 / i); 427193640Sariff Popts["rolloff"] = rolloff_round(i); 428193640Sariff return (0); 429193640Sariff } else { 430193640Sariff attn = 1.0 * a[i++]; 431193640Sariff Popts["beta"] = kaiserAttn2Beta(attn); 432193640Sariff } 433193640Sariff i = tap_round(a[i]); 434193640Sariff Popts["nmult"] = i; 435193640Sariff if (i > 7 && i < 28) 436193640Sariff i = 27; 437193640Sariff i = kaiserRolloff(i, attn); 438193640Sariff Popts["rolloff"] = rolloff_round(i); 439193640Sariff 440193640Sariff return (0); 441193640Sariff } 442193640Sariff 443193640Sariff if (!(alen == 3 || alen == 4)) 444193640Sariff return (-1); 445193640Sariff 446193640Sariff i = 2; 447193640Sariff 448193640Sariff if (a[1] == "kaiser") { 449193640Sariff if (alen > 2) 450193640Sariff Popts["beta"] = 1.0 * a[i++]; 451193640Sariff else 452193640Sariff Popts["beta"] = Z_KAISER_BETA_DEFAULT; 453193640Sariff } else if (Z_WINDOWS[a[1]] < Z_WINDOW_KAISER) 454193640Sariff Popts["beta"] = Z_WINDOWS[a[1]]; 455193640Sariff else if (1.0 * a[1] < Z_WINDOW_KAISER) 456193640Sariff return (-1); 457193640Sariff else 458193640Sariff Popts["beta"] = kaiserAttn2Beta(1.0 * a[1]); 459193640Sariff Popts["nmult"] = tap_round(a[i++]); 460193640Sariff if (a[1] == "kaiser" && alen == 3) 461193640Sariff i = kaiserRolloff(Popts["nmult"], \ 462193640Sariff kaiserBeta2Attn(Popts["beta"])); 463193640Sariff else 464193640Sariff i = 1.0 * a[i]; 465193640Sariff Popts["rolloff"] = rolloff_round(i); 466193640Sariff 467193640Sariff return (0); 468193640Sariff} 469193640Sariff 470193640Sarifffunction genscale(bit, s1, s2, scale) 471193640Sariff{ 472195378Sariff if ((bit + Z_COEFF_SHIFT) > Z_ACCUMULATOR_BIT) 473195378Sariff s1 = Z_COEFF_SHIFT - \ 474195378Sariff (32 - (Z_ACCUMULATOR_BIT - Z_COEFF_SHIFT)); 475195378Sariff else 476195378Sariff s1 = Z_COEFF_SHIFT - (32 - bit); 477195378Sariff 478193640Sariff s2 = Z_SHIFT + (32 - bit); 479193640Sariff 480193640Sariff if (s1 == 0) 481193640Sariff scale = "v"; 482193640Sariff else if (s1 < 0) 483193640Sariff scale = sprintf("(v) << %d", abs(s1)); 484193640Sariff else 485193640Sariff scale = sprintf("(v) >> %d", s1); 486193640Sariff 487193640Sariff scale = sprintf("(%s) * Z_SCALE_CAST(s)", scale); 488193640Sariff 489193640Sariff if (s2 != 0) 490193640Sariff scale = sprintf("(%s) >> %d", scale, s2); 491193640Sariff 492193640Sariff printf("#define Z_SCALE_%d(v, s)\t%s(%s)\n", \ 493193640Sariff bit, (bit < 10) ? "\t" : "", scale); 494193640Sariff} 495193640Sariff 496193640Sarifffunction genlerp(bit, use64, lerp) 497193640Sariff{ 498193640Sariff if ((bit + Z_LINEAR_SHIFT) <= 32) { 499193640Sariff lerp = sprintf("(((y) - (x)) * (z)) >> %d", Z_LINEAR_SHIFT); 500193640Sariff } else if (use64 != 0) { 501193640Sariff if ((bit + Z_LINEAR_SHIFT) <= 64) { 502193640Sariff lerp = sprintf( \ 503193640Sariff "(((int64_t)(y) - (x)) * (z)) " \ 504193640Sariff ">> %d", \ 505193640Sariff Z_LINEAR_SHIFT); 506193640Sariff } else { 507193640Sariff lerp = sprintf( \ 508193640Sariff "((int64_t)((y) >> %d) - ((x) >> %d)) * ", \ 509193640Sariff "(z)" \ 510193640Sariff bit + Z_LINEAR_SHIFT - 64, \ 511193640Sariff bit + Z_LINEAR_SHIFT - 64); 512193640Sariff if ((64 - bit) != 0) 513193640Sariff lerp = sprintf("(%s) >> %d", lerp, 64 - bit); 514193640Sariff } 515193640Sariff } else { 516193640Sariff lerp = sprintf( \ 517193640Sariff "(((y) >> %d) - ((x) >> %d)) * (z)", \ 518193640Sariff bit + Z_LINEAR_SHIFT - 32, \ 519193640Sariff bit + Z_LINEAR_SHIFT - 32); 520193640Sariff if ((32 - bit) != 0) 521193640Sariff lerp = sprintf("(%s) >> %d", lerp, 32 - bit); 522193640Sariff } 523193640Sariff 524193640Sariff printf("#define Z_LINEAR_INTERPOLATE_%d(z, x, y)" \ 525193640Sariff "\t\t\t\t%s\\\n\t((x) + (%s))\n", \ 526193640Sariff bit, (bit < 10) ? "\t" : "", lerp); 527193640Sariff} 528193640Sariff 529194233Sarifffunction init_drift(drift, xdrift) 530194233Sariff{ 531194233Sariff xdrift = floor(drift); 532194233Sariff 533194233Sariff if (Z_DRIFT_SHIFT != -1) { 534194233Sariff if (xdrift != Z_DRIFT_SHIFT) 535194233Sariff printf("#error Z_DRIFT_SHIFT reinitialize!\n"); 536194233Sariff return; 537194233Sariff } 538194233Sariff 539194233Sariff # 540194233Sariff # Initialize filter oversampling factor, or in other word 541194233Sariff # Z_DRIFT_SHIFT. 542194233Sariff # 543194233Sariff if (xdrift < 0) 544194233Sariff xdrift = 1; 545194233Sariff else if (xdrift > 31) 546194233Sariff xdrift = 31; 547194233Sariff 548194233Sariff Z_DRIFT_SHIFT = xdrift; 549194233Sariff Z_DRIFT_ONE = shl(1, Z_DRIFT_SHIFT); 550194233Sariff 551194233Sariff Z_SHIFT = Z_FULL_SHIFT - Z_DRIFT_SHIFT; 552194233Sariff Z_ONE = shl(1, Z_SHIFT); 553194233Sariff Z_MASK = Z_ONE - 1; 554194233Sariff} 555194233Sariff 556195378Sarifffunction init_coeff_bit(cbit, xcbit) 557195378Sariff{ 558195378Sariff xcbit = floor(cbit); 559195378Sariff 560195378Sariff if (Z_COEFF_SHIFT != 0) { 561195378Sariff if (xcbit != Z_COEFF_SHIFT) 562195378Sariff printf("#error Z_COEFF_SHIFT reinitialize!\n"); 563195378Sariff return; 564195378Sariff } 565195378Sariff 566195378Sariff # 567195378Sariff # Initialize dynamic range of coefficients. 568195378Sariff # 569195378Sariff if (xcbit < 1) 570195378Sariff xcbit = 1; 571195378Sariff else if (xcbit > 30) 572195378Sariff xcbit = 30; 573195378Sariff 574195378Sariff Z_COEFF_SHIFT = xcbit; 575195378Sariff Z_COEFF_ONE = shl(1, Z_COEFF_SHIFT); 576195378Sariff} 577195378Sariff 578195378Sarifffunction init_accum_bit(accbit, xaccbit) 579195378Sariff{ 580195378Sariff xaccbit = floor(accbit); 581195378Sariff 582195378Sariff if (Z_ACCUMULATOR_BIT != 0) { 583195378Sariff if (xaccbit != Z_ACCUMULATOR_BIT) 584195378Sariff printf("#error Z_ACCUMULATOR_BIT reinitialize!\n"); 585195378Sariff return; 586195378Sariff } 587195378Sariff 588195378Sariff # 589195378Sariff # Initialize dynamic range of accumulator. 590195378Sariff # 591195378Sariff if (xaccbit > 64) 592195378Sariff xaccbit = 64; 593195378Sariff else if (xaccbit < 32) 594195378Sariff xaccbit = 32; 595195378Sariff 596195378Sariff Z_ACCUMULATOR_BIT = xaccbit; 597195378Sariff} 598195378Sariff 599195378Sarifffunction init_coeff_interpolator(interp) 600195378Sariff{ 601195378Sariff # 602195378Sariff # Validate interpolator type. 603195378Sariff # 604195378Sariff if (interp == "ZOH" || interp == "LINEAR" || \ 605195378Sariff interp == "QUADRATIC" || interp == "HERMITE" || \ 606195378Sariff interp == "BSPLINE" || interp == "OPT32X" || \ 607195378Sariff interp == "OPT16X" || interp == "OPT8X" || \ 608195378Sariff interp == "OPT4X" || interp == "OPT2X") 609195378Sariff Z_COEFF_INTERPOLATOR = interp; 610195378Sariff} 611195378Sariff 612193640SariffBEGIN { 613193640Sariff I0_EPSILON = 1e-21; 614193640Sariff M_PI = atan2(0.0, -1.0); 615193640Sariff 616193640Sariff INT32_MAX = 1 + ((shl(1, 30) - 1) * 2); 617193640Sariff INT32_MIN = -1 - INT32_MAX; 618193640Sariff 619193640Sariff Z_COEFF_OFFSET = 5; 620193640Sariff 621195378Sariff Z_ACCUMULATOR_BIT_DEFAULT = 58; 622195378Sariff Z_ACCUMULATOR_BIT = 0; 623195378Sariff 624193640Sariff Z_FULL_SHIFT = 30; 625193640Sariff Z_FULL_ONE = shl(1, Z_FULL_SHIFT); 626193640Sariff 627195378Sariff Z_COEFF_SHIFT_DEFAULT = 30; 628195378Sariff Z_COEFF_SHIFT = 0; 629195378Sariff Z_COEFF_ONE = 0; 630193640Sariff 631195378Sariff Z_COEFF_INTERPOLATOR = 0; 632195378Sariff 633193640Sariff Z_INTERP_COEFF_SHIFT = 24; 634193640Sariff Z_INTERP_COEFF_ONE = shl(1, Z_INTERP_COEFF_SHIFT); 635193640Sariff 636193640Sariff Z_LINEAR_FULL_SHIFT = Z_FULL_SHIFT; 637193640Sariff Z_LINEAR_FULL_ONE = shl(1, Z_LINEAR_FULL_SHIFT); 638193640Sariff Z_LINEAR_SHIFT = 8; 639193640Sariff Z_LINEAR_UNSHIFT = Z_LINEAR_FULL_SHIFT - Z_LINEAR_SHIFT; 640193640Sariff Z_LINEAR_ONE = shl(1, Z_LINEAR_SHIFT) 641193640Sariff 642194233Sariff Z_DRIFT_SHIFT_DEFAULT = 5; 643194233Sariff Z_DRIFT_SHIFT = -1; 644193640Sariff # meehhhh... let it overflow... 645193640Sariff #Z_SCALE_SHIFT = 31; 646193640Sariff #Z_SCALE_ONE = shl(1, Z_SCALE_SHIFT); 647193640Sariff 648193640Sariff Z_WINDOW_KAISER = 0.0; 649193640Sariff Z_WINDOW_BLACKMAN_NUTTALL = -1.0; 650193640Sariff Z_WINDOW_NUTTALL = -2.0; 651193640Sariff Z_WINDOW_BLACKMAN_HARRIS = -3.0; 652193640Sariff Z_WINDOW_BLACKMAN = -4.0; 653193640Sariff Z_WINDOW_HAMMING = -5.0; 654193640Sariff Z_WINDOW_HANN = -6.0; 655193640Sariff 656193640Sariff Z_WINDOWS["blackman_nuttall"] = Z_WINDOW_BLACKMAN_NUTTALL; 657193640Sariff Z_WINDOWS["nuttall"] = Z_WINDOW_NUTTALL; 658193640Sariff Z_WINDOWS["blackman_harris"] = Z_WINDOW_BLACKMAN_HARRIS; 659193640Sariff Z_WINDOWS["blackman"] = Z_WINDOW_BLACKMAN; 660193640Sariff Z_WINDOWS["hamming"] = Z_WINDOW_HAMMING; 661193640Sariff Z_WINDOWS["hann"] = Z_WINDOW_HANN; 662193640Sariff 663193640Sariff Z_KAISER_2_BLACKMAN_BETA = 8.568611; 664193640Sariff Z_KAISER_2_BLACKMAN_NUTTALL_BETA = 11.98; 665193640Sariff 666193640Sariff Z_KAISER_ATTN_DEFAULT = 100; 667193640Sariff Z_KAISER_BETA_DEFAULT = kaiserAttn2Beta(Z_KAISER_ATTN_DEFAULT); 668193640Sariff 669193640Sariff Z_KAISER_EPSILON = 1e-21; 670193640Sariff 671193640Sariff # 672193640Sariff # This is practically a joke. 673193640Sariff # 674193640Sariff Z_NYQUIST_HOVER = 0.0; 675193640Sariff 676193640Sariff smallest = 10.000000; 677193640Sariff largest = 0.000010; 678193640Sariff largest_interp = 0; 679193640Sariff 680193640Sariff if (ARGC < 2) { 681193640Sariff ARGC = 1; 682193640Sariff ARGV[ARGC++] = "100:8:0.85"; 683195283Sariff ARGV[ARGC++] = "100:36:0.92"; 684193640Sariff ARGV[ARGC++] = "100:164:0.97"; 685193640Sariff #ARGV[ARGC++] = "100:8"; 686193640Sariff #ARGV[ARGC++] = "100:16"; 687193640Sariff #ARGV[ARGC++] = "100:32:0.7929"; 688193640Sariff #ARGV[ARGC++] = "100:64:0.8990"; 689193640Sariff #ARGV[ARGC++] = "100:128:0.9499"; 690193640Sariff } 691193640Sariff 692193640Sariff printf("#ifndef _FEEDER_RATE_GEN_H_\n"); 693193640Sariff printf("#define _FEEDER_RATE_GEN_H_\n\n"); 694193640Sariff printf("/*\n"); 695193640Sariff printf(" * Generated using feeder_rate_mkfilter.awk, heaven, wind and awesome.\n"); 696193640Sariff printf(" *\n"); 697193640Sariff printf(" * DO NOT EDIT!\n"); 698193640Sariff printf(" */\n\n"); 699193640Sariff printf("#define FEEDER_RATE_PRESETS\t\""); 700193640Sariff for (i = 1; i < ARGC; i++) 701193640Sariff printf("%s%s", (i == 1) ? "" : " ", ARGV[i]); 702193640Sariff printf("\"\n\n"); 703193640Sariff imp["quality"] = 2; 704193640Sariff for (i = 1; i < ARGC; i++) { 705193640Sariff if (filter_parse(ARGV[i]) == 0) { 706193640Sariff beta = Popts["beta"]; 707193640Sariff nmult = Popts["nmult"]; 708193640Sariff rolloff = Popts["rolloff"]; 709194233Sariff if (Z_DRIFT_SHIFT == -1) 710194233Sariff init_drift(Z_DRIFT_SHIFT_DEFAULT); 711195378Sariff if (Z_COEFF_SHIFT == 0) 712195378Sariff init_coeff_bit(Z_COEFF_SHIFT_DEFAULT); 713195378Sariff if (Z_ACCUMULATOR_BIT == 0) 714195378Sariff init_accum_bit(Z_ACCUMULATOR_BIT_DEFAULT); 715193640Sariff ztab[imp["quality"] - 2] = \ 716193640Sariff mkfilter(imp, nmult, rolloff, beta, Z_DRIFT_ONE); 717193640Sariff imp["quality"]++; 718193640Sariff } 719193640Sariff } 720193640Sariff 721193640Sariff printf("\n"); 722193640Sariff # 723193640Sariff # XXX 724193640Sariff # 725193640Sariff #if (length(ztab) > 0) { 726193640Sariff # j = 0; 727193640Sariff # for (i = 0; i < length(ztab); i++) { 728193640Sariff # if (ztab[i] > j) 729193640Sariff # j = ztab[i]; 730193640Sariff # } 731193640Sariff # printf("static int32_t z_coeff_zero[%d] = {", j); 732193640Sariff # for (i = 0; i < j; i++) { 733193640Sariff # if ((i % 19) == 0) 734193640Sariff # printf("\n"); 735193640Sariff # printf(" 0,"); 736193640Sariff # } 737193640Sariff # printf("\n};\n\n"); 738193640Sariff #} 739193640Sariff # 740193640Sariff # XXX 741193640Sariff # 742193640Sariff printf("static const struct {\n"); 743193640Sariff printf("\tint32_t len;\n"); 744193640Sariff printf("\tint32_t *coeff;\n"); 745193640Sariff printf("\tint32_t *dcoeff;\n"); 746193640Sariff printf("} z_coeff_tab[] = {\n"); 747193640Sariff if (length(ztab) > 0) { 748193640Sariff j = 0; 749193640Sariff for (i = 0; i < length(ztab); i++) { 750193640Sariff if (ztab[i] > j) 751193640Sariff j = ztab[i]; 752193640Sariff } 753193640Sariff j = length(sprintf("%d", j)); 754193640Sariff lfmt = sprintf("%%%dd", j); 755193640Sariff j = length(sprintf("z_coeff_q%d", length(ztab) + 1)); 756193640Sariff zcfmt = sprintf("%%-%ds", j); 757193640Sariff zdcfmt = sprintf("%%-%ds", j + 1); 758193640Sariff 759193640Sariff for (i = 0; i < length(ztab); i++) { 760193640Sariff l = sprintf(lfmt, ztab[i]); 761193640Sariff zc = sprintf("z_coeff_q%d", i + 2); 762193640Sariff zc = sprintf(zcfmt, zc); 763193640Sariff zdc = sprintf("z_dcoeff_q%d", i + 2); 764193640Sariff zdc = sprintf(zdcfmt, zdc); 765193640Sariff printf("\t{ %s, %s, %s },\n", l, zc, zdc); 766193640Sariff } 767193640Sariff } else 768193640Sariff printf("\t{ 0, NULL, NULL }\n"); 769193640Sariff printf("};\n\n"); 770193640Sariff 771193640Sariff #Z_UNSHIFT = 0; 772193640Sariff #v = shr(Z_ONE - 1, Z_UNSHIFT) * abs(largest_interp); 773193640Sariff #while (v < 0 || v > INT32_MAX) { 774193640Sariff # Z_UNSHIFT += 1; 775193640Sariff # v = shr(Z_ONE - 1, Z_UNSHIFT) * abs(largest_interp); 776193640Sariff #} 777193640Sariff v = ((Z_ONE - 1) * abs(largest_interp)) / INT32_MAX; 778193640Sariff Z_UNSHIFT = ceil(log(v) / log(2.0)); 779193640Sariff Z_INTERP_SHIFT = Z_SHIFT - Z_UNSHIFT + Z_INTERP_COEFF_SHIFT; 780193640Sariff 781193640Sariff Z_INTERP_UNSHIFT = (Z_SHIFT - Z_UNSHIFT) + Z_INTERP_COEFF_SHIFT \ 782193640Sariff - Z_COEFF_SHIFT; 783193640Sariff 784193640Sariff printf("#define Z_COEFF_TAB_SIZE\t\t\t\t\t\t\\\n"); 785193640Sariff printf("\t((int32_t)(sizeof(z_coeff_tab) /"); 786193640Sariff printf(" sizeof(z_coeff_tab[0])))\n\n"); 787193640Sariff printf("#define Z_COEFF_OFFSET\t\t%d\n\n", Z_COEFF_OFFSET); 788193640Sariff printf("#define Z_RSHIFT(x, y)\t\t(((x) + " \ 789193640Sariff "(1 << ((y) - 1))) >> (y))\n"); 790193640Sariff printf("#define Z_RSHIFT_L(x, y)\t(((x) + " \ 791193640Sariff "(1LL << ((y) - 1))) >> (y))\n\n"); 792193640Sariff printf("#define Z_FULL_SHIFT\t\t%d\n", Z_FULL_SHIFT); 793193640Sariff printf("#define Z_FULL_ONE\t\t0x%08x%s\n", Z_FULL_ONE, \ 794193640Sariff (Z_FULL_ONE > INT32_MAX) ? "U" : ""); 795193640Sariff printf("\n"); 796193640Sariff printf("#define Z_DRIFT_SHIFT\t\t%d\n", Z_DRIFT_SHIFT); 797193640Sariff #printf("#define Z_DRIFT_ONE\t\t0x%08x\n", Z_DRIFT_ONE); 798193640Sariff printf("\n"); 799193640Sariff printf("#define Z_SHIFT\t\t\t%d\n", Z_SHIFT); 800193640Sariff printf("#define Z_ONE\t\t\t0x%08x\n", Z_ONE); 801193640Sariff printf("#define Z_MASK\t\t\t0x%08x\n", Z_MASK); 802193640Sariff printf("\n"); 803193640Sariff printf("#define Z_COEFF_SHIFT\t\t%d\n", Z_COEFF_SHIFT); 804193640Sariff zinterphp = "(z) * (d)"; 805193640Sariff zinterpunshift = Z_SHIFT + Z_INTERP_COEFF_SHIFT - Z_COEFF_SHIFT; 806193640Sariff if (zinterpunshift > 0) { 807193640Sariff v = (Z_ONE - 1) * abs(largest_interp); 808193640Sariff if (v < INT32_MIN || v > INT32_MAX) 809193640Sariff zinterphp = sprintf("(int64_t)%s", zinterphp); 810193640Sariff zinterphp = sprintf("(%s) >> %d", zinterphp, zinterpunshift); 811193640Sariff } else if (zinterpunshift < 0) 812193640Sariff zinterphp = sprintf("(%s) << %d", zinterphp, \ 813193640Sariff abs(zinterpunshift)); 814193640Sariff if (Z_UNSHIFT == 0) 815193640Sariff zinterp = "z"; 816193640Sariff else 817193640Sariff zinterp = sprintf("(z) >> %d", Z_UNSHIFT); 818193640Sariff zinterp = sprintf("(%s) * (d)", zinterp); 819193640Sariff if (Z_INTERP_UNSHIFT < 0) 820193640Sariff zinterp = sprintf("(%s) << %d", zinterp, \ 821193640Sariff abs(Z_INTERP_UNSHIFT)); 822193640Sariff else if (Z_INTERP_UNSHIFT > 0) 823193640Sariff zinterp = sprintf("(%s) >> %d", zinterp, Z_INTERP_UNSHIFT); 824193640Sariff if (zinterphp != zinterp) { 825193640Sariff printf("\n#ifdef SND_FEEDER_RATE_HP\n"); 826193640Sariff printf("#define Z_COEFF_INTERPOLATE(z, c, d)" \ 827193640Sariff "\t\t\t\t\t\\\n\t((c) + (%s))\n", zinterphp); 828193640Sariff printf("#else\n"); 829193640Sariff printf("#define Z_COEFF_INTERPOLATE(z, c, d)" \ 830193640Sariff "\t\t\t\t\t\\\n\t((c) + (%s))\n", zinterp); 831193640Sariff printf("#endif\n"); 832193640Sariff } else 833193640Sariff printf("#define Z_COEFF_INTERPOLATE(z, c, d)" \ 834193640Sariff "\t\t\t\t\t\\\n\t((c) + (%s))\n", zinterp); 835193640Sariff #printf("\n"); 836193640Sariff #printf("#define Z_SCALE_SHIFT\t\t%d\n", Z_SCALE_SHIFT); 837193640Sariff #printf("#define Z_SCALE_ONE\t\t0x%08x%s\n", Z_SCALE_ONE, \ 838193640Sariff # (Z_SCALE_ONE > INT32_MAX) ? "U" : ""); 839193640Sariff printf("\n"); 840193640Sariff printf("#define Z_SCALE_CAST(s)\t\t((uint32_t)(s))\n"); 841193640Sariff genscale(8); 842193640Sariff genscale(16); 843193640Sariff genscale(24); 844193640Sariff genscale(32); 845193640Sariff printf("\n"); 846195378Sariff printf("#define Z_ACCUMULATOR_BIT\t%d\n\n", Z_ACCUMULATOR_BIT) 847195378Sariff for (i = 8; i <= 32; i += 8) { 848195378Sariff gbit = ((i + Z_COEFF_SHIFT) > Z_ACCUMULATOR_BIT) ? \ 849195378Sariff (i - (Z_ACCUMULATOR_BIT - Z_COEFF_SHIFT)) : 0; 850195378Sariff printf("#define Z_GUARD_BIT_%d\t\t%d\n", i, gbit); 851195378Sariff if (gbit == 0) 852195378Sariff printf("#define Z_NORM_%d(v)\t\t(v)\n\n", i); 853195378Sariff else 854195378Sariff printf("#define Z_NORM_%d(v)\t\t" \ 855195378Sariff "((v) >> Z_GUARD_BIT_%d)\n\n", i, i); 856195378Sariff } 857195378Sariff printf("\n"); 858193640Sariff printf("#define Z_LINEAR_FULL_ONE\t0x%08xU\n", Z_LINEAR_FULL_ONE); 859193640Sariff printf("#define Z_LINEAR_SHIFT\t\t%d\n", Z_LINEAR_SHIFT); 860193640Sariff printf("#define Z_LINEAR_UNSHIFT\t%d\n", Z_LINEAR_UNSHIFT); 861193640Sariff printf("#define Z_LINEAR_ONE\t\t0x%08x\n", Z_LINEAR_ONE); 862193640Sariff printf("\n"); 863193640Sariff printf("#ifdef SND_PCM_64\n"); 864193640Sariff genlerp(8, 1); 865193640Sariff genlerp(16, 1); 866193640Sariff genlerp(24, 1); 867193640Sariff genlerp(32, 1); 868193640Sariff printf("#else\t/* !SND_PCM_64 */\n"); 869193640Sariff genlerp(8, 0); 870193640Sariff genlerp(16, 0); 871193640Sariff genlerp(24, 0); 872193640Sariff genlerp(32, 0); 873193640Sariff printf("#endif\t/* SND_PCM_64 */\n"); 874193640Sariff printf("\n"); 875193640Sariff printf("#define Z_QUALITY_ZOH\t\t0\n"); 876193640Sariff printf("#define Z_QUALITY_LINEAR\t1\n"); 877193640Sariff printf("#define Z_QUALITY_SINC\t\t%d\n", \ 878193640Sariff floor((length(ztab) - 1) / 2) + 2); 879193640Sariff printf("\n"); 880193640Sariff printf("#define Z_QUALITY_MIN\t\t0\n"); 881193640Sariff printf("#define Z_QUALITY_MAX\t\t%d\n", length(ztab) + 1); 882195378Sariff if (Z_COEFF_INTERPOLATOR != 0) 883195378Sariff printf("\n#define Z_COEFF_INTERP_%s\t1\n", \ 884195378Sariff Z_COEFF_INTERPOLATOR); 885193640Sariff printf("\n/*\n * smallest: %.32f\n * largest: %.32f\n *\n", \ 886193640Sariff smallest, largest); 887193640Sariff printf(" * z_unshift=%d, z_interp_shift=%d\n *\n", \ 888193640Sariff Z_UNSHIFT, Z_INTERP_SHIFT); 889193640Sariff v = shr(Z_ONE - 1, Z_UNSHIFT) * abs(largest_interp); 890193640Sariff printf(" * largest interpolation multiplication: %d\n */\n", v); 891193640Sariff if (v < INT32_MIN || v > INT32_MAX) { 892193640Sariff printf("\n#ifndef SND_FEEDER_RATE_HP\n"); 893193640Sariff printf("#error interpolation overflow, please reduce" \ 894193640Sariff " Z_INTERP_SHIFT\n"); 895193640Sariff printf("#endif\n"); 896193640Sariff } 897193640Sariff 898193640Sariff printf("\n#endif /* !_FEEDER_RATE_GEN_H_ */\n"); 899193640Sariff} 900