Alc/nfcfilter.c

   1
   2 #include "config.h"
   3
   4 #include "nfcfilter.h"
   5
   6 #include "alu.h"
   7
   8
   9 /* Near-field control filters are the basis for handling the near-field effect.
  10  * The near-field effect is a bass-boost present in the directional components
  11  * of a recorded signal, created as a result of the wavefront curvature (itself
  12  * a function of sound distance). Proper reproduction dictates this be
  13  * compensated for using a bass-cut given the playback speaker distance, to
  14  * avoid excessive bass in the playback.
  15  *
  16  * For real-time rendered audio, emulating the near-field effect based on the
  17  * sound source's distance, and subsequently compensating for it at output
  18  * based on the speaker distances, can create a more realistic perception of
  19  * sound distance beyond a simple 1/r attenuation.
  20  *
  21  * These filters do just that. Each one applies a low-shelf filter, created as
  22  * the combination of a bass-boost for a given sound source distance (near-
  23  * field emulation) along with a bass-cut for a given control/speaker distance
  24  * (near-field compensation).
  25  *
  26  * Note that it is necessary to apply a cut along with the boost, since the
  27  * boost alone is unstable in higher-order ambisonics as it causes an infinite
  28  * DC gain (even first-order ambisonics requires there to be no DC offset for
  29  * the boost to work). Consequently, ambisonics requires a control parameter to
  30  * be used to avoid an unstable boost-only filter. NFC-HOA defines this control
  31  * as a reference delay, calculated with:
  32  *
  33  * reference_delay = control_distance / speed_of_sound
  34  *
  35  * This means w0 (for input) or w1 (for output) should be set to:
  36  *
  37  * wN = 1 / (reference_delay * sample_rate)
  38  *
  39  * when dealing with NFC-HOA content. For FOA input content, which does not
  40  * specify a reference_delay variable, w0 should be set to 0 to apply only
  41  * near-field compensation for output. It's important that w1 be a finite,
  42  * positive, non-0 value or else the bass-boost will become unstable again.
  43  * Also, w0 should not be too large compared to w1, to avoid excessively loud
  44  * low frequencies.
  45  */
  46
  47 static const float B[4][3] = {
  48     {    0.0f                             },
  49     {    1.0f                             },
  50     {    3.0f,     3.0f                   },
  51     { 3.6778f,  6.4595f, 2.3222f          },
  52   /*{ 4.2076f, 11.4877f, 5.7924f, 9.1401f }*/
  53 };
  54
  55 void NfcFilterCreate1(NfcFilter *nfc, const float w0, const float w1)
  56 {
  57     float b_00, g_0;
  58     float r;
  59
  60     memset(nfc, 0, sizeof(*nfc));
  61
  62     nfc->g = 1.0f;
  63     nfc->coeffs[0] = 1.0f;
  64
  65     /* Calculate bass-boost coefficients. */
  66     r = 0.5f * w0;
  67     b_00 = B[1][0] * r;
  68     g_0 = 1.0f + b_00;
  69
  70     nfc->coeffs[0] *= g_0;
  71     nfc->coeffs[1] = (2.0f * b_00) / g_0;
  72
  73     /* Calculate bass-cut coefficients. */
  74     r = 0.5f * w1;
  75     b_00 = B[1][0] * r;
  76     g_0 = 1.0f + b_00;
  77
  78     nfc->g /= g_0;
  79     nfc->coeffs[0] /= g_0;
  80     nfc->coeffs[1+1] = (2.0f * b_00) / g_0;
  81 }
  82
  83 void NfcFilterAdjust1(NfcFilter *nfc, const float w0)
  84 {
  85     float b_00, g_0;
  86     float r;
  87
  88     r = 0.5f * w0;
  89     b_00 = B[1][0] * r;
  90     g_0 = 1.0f + b_00;
  91
  92     nfc->coeffs[0] = nfc->g * g_0;
  93     nfc->coeffs[1] = (2.0f * b_00) / g_0;
  94 }
  95
  96 void NfcFilterUpdate1(NfcFilter *nfc, ALfloat *restrict dst, const float *restrict src, const int count)
  97 {
  98     const float b0 = nfc->coeffs[0];
  99     const float a0 = nfc->coeffs[1];
 100     const float a1 = nfc->coeffs[2];
 101     float z1 = nfc->history[0];
 102     int i;
 103
 104     for(i = 0;i < count;i++)
 105     {
 106         float out = src[i] * b0;
 107         float y;
 108
 109         y = out - (a1*z1);
 110         out = y + (a0*z1);
 111         z1 += y;
 112
 113         dst[i] = out;
 114     }
 115     nfc->history[0] = z1;
 116 }
 117
 118
 119 void NfcFilterCreate2(NfcFilter *nfc, const float w0, const float w1)
 120 {
 121     float b_10, b_11, g_1;
 122     float r;
 123
 124     memset(nfc, 0, sizeof(*nfc));
 125
 126     nfc->g = 1.0f;
 127     nfc->coeffs[0] = 1.0f;
 128
 129     /* Calculate bass-boost coefficients. */
 130     r = 0.5f * w0;
 131     b_10 = B[2][0] * r;
 132     b_11 = B[2][1] * r * r;
 133     g_1 = 1.0f + b_10 + b_11;
 134
 135     nfc->coeffs[0] *= g_1;
 136     nfc->coeffs[1] = ((2.0f * b_10) + (4.0f * b_11)) / g_1;
 137     nfc->coeffs[2] = (4.0f * b_11) / g_1;
 138
 139     /* Calculate bass-cut coefficients. */
 140     r = 0.5f * w1;
 141     b_10 = B[2][0] * r;
 142     b_11 = B[2][1] * r * r;
 143     g_1 = 1.0f + b_10 + b_11;
 144
 145     nfc->g /= g_1;
 146     nfc->coeffs[0] /= g_1;
 147     nfc->coeffs[2+1] = ((2.0f * b_10) + (4.0f * b_11)) / g_1;
 148     nfc->coeffs[2+2] = (4.0f * b_11) / g_1;
 149 }
 150
 151 void NfcFilterAdjust2(NfcFilter *nfc, const float w0)
 152 {
 153     float b_10, b_11, g_1;
 154     float r;
 155
 156     r = 0.5f * w0;
 157     b_10 = B[2][0] * r;
 158     b_11 = B[2][1] * r * r;
 159     g_1 = 1.0f + b_10 + b_11;
 160
 161     nfc->coeffs[0] = nfc->g * g_1;
 162     nfc->coeffs[1] = ((2.0f * b_10) + (4.0f * b_11)) / g_1;
 163     nfc->coeffs[2] = (4.0f * b_11) / g_1;
 164 }
 165
 166 void NfcFilterUpdate2(NfcFilter *nfc, ALfloat *restrict dst, const float *restrict src, const int count)
 167 {
 168     const float b0 = nfc->coeffs[0];
 169     const float a00 = nfc->coeffs[1];
 170     const float a01 = nfc->coeffs[2];
 171     const float a10 = nfc->coeffs[3];
 172     const float a11 = nfc->coeffs[4];
 173     float z1 = nfc->history[0];
 174     float z2 = nfc->history[1];
 175     int i;
 176
 177     for(i = 0;i < count;i++)
 178     {
 179         float out = src[i] * b0;
 180         float y;
 181
 182         y = out - (a10*z1) - (a11*z2);
 183         out = y + (a00*z1) + (a01*z2);
 184         z2 += z1;
 185         z1 += y;
 186
 187         dst[i] = out;
 188     }
 189     nfc->history[0] = z1;
 190     nfc->history[1] = z2;
 191 }
 192
 193
 194 void NfcFilterCreate3(NfcFilter *nfc, const float w0, const float w1)
 195 {
 196     float b_10, b_11, g_1;
 197     float b_00, g_0;
 198     float r;
 199
 200     memset(nfc, 0, sizeof(*nfc));
 201
 202     nfc->g = 1.0f;
 203     nfc->coeffs[0] = 1.0f;
 204
 205     /* Calculate bass-boost coefficients. */
 206     r = 0.5f * w0;
 207     b_10 = B[3][0] * r;
 208     b_11 = B[3][1] * r * r;
 209     g_1 = 1.0f + b_10 + b_11;
 210
 211     nfc->coeffs[0] *= g_1;
 212     nfc->coeffs[1] = ((2.0f * b_10) + (4.0f * b_11)) / g_1;
 213     nfc->coeffs[2] = (4.0f * b_11) / g_1;
 214
 215     b_00 = B[3][2] * r;
 216     g_0 = 1.0f + b_00;
 217
 218     nfc->coeffs[0] *= g_0;
 219     nfc->coeffs[2+1] = (2.0f * b_00) / g_0;
 220
 221     /* Calculate bass-cut coefficients. */
 222     r = 0.5f * w1;
 223     b_10 = B[3][0] * r;
 224     b_11 = B[3][1] * r * r;
 225     g_1 = 1.0f + b_10 + b_11;
 226
 227     nfc->g /= g_1;
 228     nfc->coeffs[0] /= g_1;
 229     nfc->coeffs[3+1] = ((2.0f * b_10) + (4.0f * b_11)) / g_1;
 230     nfc->coeffs[3+2] = (4.0f * b_11) / g_1;
 231
 232     b_00 = B[3][2] * r;
 233     g_0 = 1.0f + b_00;
 234
 235     nfc->g /= g_0;
 236     nfc->coeffs[0] /= g_0;
 237     nfc->coeffs[3+2+1] = (2.0f * b_00) / g_0;
 238 }
 239
 240 void NfcFilterAdjust3(NfcFilter *nfc, const float w0)
 241 {
 242     float b_10, b_11, g_1;
 243     float b_00, g_0;
 244     float r;
 245
 246     r = 0.5f * w0;
 247     b_10 = B[3][0] * r;
 248     b_11 = B[3][1] * r * r;
 249     g_1 = 1.0f + b_10 + b_11;
 250
 251     nfc->coeffs[0] = nfc->g * g_1;
 252     nfc->coeffs[1] = ((2.0f * b_10) + (4.0f * b_11)) / g_1;
 253     nfc->coeffs[2] = (4.0f * b_11) / g_1;
 254
 255     b_00 = B[3][2] * r;
 256     g_0 = 1.0f + b_00;
 257
 258     nfc->coeffs[0] *= g_0;
 259     nfc->coeffs[2+1] = (2.0f * b_00) / g_0;
 260 }
 261
 262 void NfcFilterUpdate3(NfcFilter *nfc, ALfloat *restrict dst, const float *restrict src, const int count)
 263 {
 264     const float b0 = nfc->coeffs[0];
 265     const float a00 = nfc->coeffs[1];
 266     const float a01 = nfc->coeffs[2];
 267     const float a02 = nfc->coeffs[3];
 268     const float a10 = nfc->coeffs[4];
 269     const float a11 = nfc->coeffs[5];
 270     const float a12 = nfc->coeffs[6];
 271     float z1 = nfc->history[0];
 272     float z2 = nfc->history[1];
 273     float z3 = nfc->history[2];
 274     int i;
 275
 276     for(i = 0;i < count;i++)
 277     {
 278         float out = src[i] * b0;
 279         float y;
 280
 281         y = out - (a10*z1) - (a11*z2);
 282         out = y + (a00*z1) + (a01*z2);
 283         z2 += z1;
 284         z1 += y;
 285
 286         y = out - (a12*z3);
 287         out = y + (a02*z3);
 288         z3 += y;
 289
 290         dst[i] = out;
 291     }
 292     nfc->history[0] = z1;
 293     nfc->history[1] = z2;
 294     nfc->history[2] = z3;
 295 }
 296
 297
 298 #if 0 /* Original methods the above are derived from. */
 299 static void NfcFilterCreate(NfcFilter *nfc, const ALsizei order, const float src_dist, const float ctl_dist, const float rate)
 300 {
 301     static const float B[4][5] = {
 302         {                                     },
 303         {    1.0f                             },
 304         {    3.0f,     3.0f                   },
 305         { 3.6778f,  6.4595f, 2.3222f          },
 306         { 4.2076f, 11.4877f, 5.7924f, 9.1401f }
 307     };
 308     float w0 = SPEEDOFSOUNDMETRESPERSEC / (src_dist * rate);
 309     float w1 = SPEEDOFSOUNDMETRESPERSEC / (ctl_dist * rate);
 310     ALsizei i;
 311     float r;
 312
 313     nfc->g = 1.0f;
 314     nfc->coeffs[0] = 1.0f;
 315
 316     /* NOTE: Slight adjustment from the literature to raise the center
 317      * frequency a bit (0.5 -> 1.0).
 318      */
 319     r = 1.0f * w0;
 320     for(i = 0; i < (order-1);i += 2)
 321     {
 322         float b_10 = B[order][i  ] * r;
 323         float b_11 = B[order][i+1] * r * r;
 324         float g_1 = 1.0f + b_10 + b_11;
 325
 326         nfc->b[i] = b_10;
 327         nfc->b[i + 1] = b_11;
 328         nfc->coeffs[0] *= g_1;
 329         nfc->coeffs[i+1] = ((2.0f * b_10) + (4.0f * b_11)) / g_1;
 330         nfc->coeffs[i+2] = (4.0f * b_11) / g_1;
 331     }
 332     if(i < order)
 333     {
 334         float b_00 = B[order][i] * r;
 335         float g_0 = 1.0f + b_00;
 336
 337         nfc->b[i] = b_00;
 338         nfc->coeffs[0] *= g_0;
 339         nfc->coeffs[i+1] = (2.0f * b_00) / g_0;
 340     }
 341
 342     r = 1.0f * w1;
 343     for(i = 0;i < (order-1);i += 2)
 344     {
 345         float b_10 = B[order][i  ] * r;
 346         float b_11 = B[order][i+1] * r * r;
 347         float g_1 = 1.0f + b_10 + b_11;
 348
 349         nfc->g /= g_1;
 350         nfc->coeffs[0] /= g_1;
 351         nfc->coeffs[order+i+1] = ((2.0f * b_10) + (4.0f * b_11)) / g_1;
 352         nfc->coeffs[order+i+2] = (4.0f * b_11) / g_1;
 353     }
 354     if(i < order)
 355     {
 356         float b_00 = B[order][i] * r;
 357         float g_0 = 1.0f + b_00;
 358
 359         nfc->g /= g_0;
 360         nfc->coeffs[0] /= g_0;
 361         nfc->coeffs[order+i+1] = (2.0f * b_00) / g_0;
 362     }
 363
 364     for(i = 0; i < MAX_AMBI_ORDER; i++)
 365         nfc->history[i] = 0.0f;
 366 }
 367
 368 static void NfcFilterAdjust(NfcFilter *nfc, const float distance)
 369 {
 370     int i;
 371
 372     nfc->coeffs[0] = nfc->g;
 373
 374     for(i = 0;i < (nfc->order-1);i += 2)
 375     {
 376         float b_10 = nfc->b[i] / distance;
 377         float b_11 = nfc->b[i+1] / (distance * distance);
 378         float g_1 = 1.0f + b_10 + b_11;
 379
 380         nfc->coeffs[0] *= g_1;
 381         nfc->coeffs[i+1] = ((2.0f * b_10) + (4.0f * b_11)) / g_1;
 382         nfc->coeffs[i+2] = (4.0f * b_11) / g_1;
 383     }
 384     if(i < nfc->order)
 385     {
 386         float b_00 = nfc->b[i] / distance;
 387         float g_0 = 1.0f + b_00;
 388
 389         nfc->coeffs[0] *= g_0;
 390         nfc->coeffs[i+1] = (2.0f * b_00) / g_0;
 391     }
 392 }
 393
 394 static float NfcFilterUpdate(const float in, NfcFilter *nfc)
 395 {
 396     int i;
 397     float out = in * nfc->coeffs[0];
 398
 399     for(i = 0;i < (nfc->order-1);i += 2)
 400     {
 401         float y = out - (nfc->coeffs[nfc->order+i+1] * nfc->history[i]) -
 402                         (nfc->coeffs[nfc->order+i+2] * nfc->history[i+1]) + 1.0e-30f;
 403         out = y + (nfc->coeffs[i+1]*nfc->history[i]) + (nfc->coeffs[i+2]*nfc->history[i+1]);
 404
 405         nfc->history[i+1] += nfc->history[i];
 406         nfc->history[i] += y;
 407     }
 408     if(i < nfc->order)
 409     {
 410         float y = out - (nfc->coeffs[nfc->order+i+1] * nfc->history[i]) + 1.0e-30f;
 411
 412         out = y + (nfc->coeffs[i+1] * nfc->history[i]);
 413         nfc->history[i] += y;
 414     }
 415
 416     return out;
 417 }
 418 #endif