Alc/filters/nfc.cpp

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