Alc/bformatdec.c

   1
   2 #include "config.h"
   3
   4 #include "bformatdec.h"
   5 #include "ambdec.h"
   6 #include "alu.h"
   7
   8 #include "threads.h"
   9
  10
  11 typedef struct BandSplitter {
  12     ALfloat coeff;
  13     ALfloat lp_z1;
  14     ALfloat lp_z2;
  15     ALfloat hp_z1;
  16 } BandSplitter;
  17
  18 static void bandsplit_init(BandSplitter *splitter, ALfloat freq_mult)
  19 {
  20     ALfloat w = freq_mult * F_TAU;
  21     ALfloat cw = cosf(w);
  22     if(cw > FLT_EPSILON)
  23         splitter->coeff = (sinf(w) - 1.0f) / cw;
  24     else
  25         splitter->coeff = cw * -0.5f;
  26
  27     splitter->lp_z1 = 0.0f;
  28     splitter->lp_z2 = 0.0f;
  29     splitter->hp_z1 = 0.0f;
  30 }
  31
  32 static void bandsplit_process(BandSplitter *splitter, ALfloat *restrict hpout, ALfloat *restrict lpout,
  33                               const ALfloat *input, ALuint count)
  34 {
  35     ALfloat coeff, d, x;
  36     ALuint i;
  37
  38     coeff = splitter->coeff*0.5f + 0.5f;
  39     for(i = 0;i < count;i++)
  40     {
  41         x = input[i];
  42
  43         d = (x - splitter->lp_z1) * coeff;
  44         x = splitter->lp_z1 + d;
  45         splitter->lp_z1 = x + d;
  46
  47         d = (x - splitter->lp_z2) * coeff;
  48         x = splitter->lp_z2 + d;
  49         splitter->lp_z2 = x + d;
  50
  51         lpout[i] = x;
  52     }
  53
  54     coeff = splitter->coeff;
  55     for(i = 0;i < count;i++)
  56     {
  57         x = input[i];
  58
  59         d = x - coeff*splitter->hp_z1;
  60         x = splitter->hp_z1 + coeff*d;
  61         splitter->hp_z1 = d;
  62
  63         hpout[i] = x - lpout[i];
  64     }
  65 }
  66
  67
  68 static const ALfloat UnitScale[MAX_AMBI_COEFFS] = {
  69     1.0f, 1.0f, 1.0f, 1.0f, 1.0f, 1.0f, 1.0f, 1.0f,
  70     1.0f, 1.0f, 1.0f, 1.0f, 1.0f, 1.0f, 1.0f, 1.0f
  71 };
  72 static const ALfloat SN3D2N3DScale[MAX_AMBI_COEFFS] = {
  73     1.000000000f, /* ACN  0 (W), sqrt(1) */
  74     1.732050808f, /* ACN  1 (Y), sqrt(3) */
  75     1.732050808f, /* ACN  2 (Z), sqrt(3) */
  76     1.732050808f, /* ACN  3 (X), sqrt(3) */
  77     2.236067978f, /* ACN  4 (V), sqrt(5) */
  78     2.236067978f, /* ACN  5 (T), sqrt(5) */
  79     2.236067978f, /* ACN  6 (R), sqrt(5) */
  80     2.236067978f, /* ACN  7 (S), sqrt(5) */
  81     2.236067978f, /* ACN  8 (U), sqrt(5) */
  82     2.645751311f, /* ACN  9 (Q), sqrt(7) */
  83     2.645751311f, /* ACN 10 (O), sqrt(7) */
  84     2.645751311f, /* ACN 11 (M), sqrt(7) */
  85     2.645751311f, /* ACN 12 (K), sqrt(7) */
  86     2.645751311f, /* ACN 13 (L), sqrt(7) */
  87     2.645751311f, /* ACN 14 (N), sqrt(7) */
  88     2.645751311f, /* ACN 15 (P), sqrt(7) */
  89 };
  90 static const ALfloat FuMa2N3DScale[MAX_AMBI_COEFFS] = {
  91     1.414213562f, /* ACN  0 (W), sqrt(2) */
  92     1.732050808f, /* ACN  1 (Y), sqrt(3) */
  93     1.732050808f, /* ACN  2 (Z), sqrt(3) */
  94     1.732050808f, /* ACN  3 (X), sqrt(3) */
  95     1.936491673f, /* ACN  4 (V), sqrt(15)/2 */
  96     1.936491673f, /* ACN  5 (T), sqrt(15)/2 */
  97     2.236067978f, /* ACN  6 (R), sqrt(5) */
  98     1.936491673f, /* ACN  7 (S), sqrt(15)/2 */
  99     1.936491673f, /* ACN  8 (U), sqrt(15)/2 */
 100     2.091650066f, /* ACN  9 (Q), sqrt(35/8) */
 101     1.972026594f, /* ACN 10 (O), sqrt(35)/3 */
 102     2.231093404f, /* ACN 11 (M), sqrt(224/45) */
 103     2.645751311f, /* ACN 12 (K), sqrt(7) */
 104     2.231093404f, /* ACN 13 (L), sqrt(224/45) */
 105     1.972026594f, /* ACN 14 (N), sqrt(35)/3 */
 106     2.091650066f, /* ACN 15 (P), sqrt(35/8) */
 107 };
 108
 109
 110 static const ALfloat SquareMatrixHF[4][MAX_AMBI_COEFFS] = {
 111     { 0.353553f,  0.204094f,  0.0f,  0.204094f },
 112     { 0.353553f, -0.204094f,  0.0f,  0.204094f },
 113     { 0.353553f,  0.204094f,  0.0f, -0.204094f },
 114     { 0.353553f, -0.204094f,  0.0f, -0.204094f },
 115 };
 116 static ALfloat SquareEncoder[4][MAX_AMBI_COEFFS];
 117
 118 static const ALfloat CubeMatrixHF[8][MAX_AMBI_COEFFS] = {
 119     { 0.25f,  0.14425f,  0.14425f,  0.14425f },
 120     { 0.25f, -0.14425f,  0.14425f,  0.14425f },
 121     { 0.25f,  0.14425f,  0.14425f, -0.14425f },
 122     { 0.25f, -0.14425f,  0.14425f, -0.14425f },
 123     { 0.25f,  0.14425f, -0.14425f,  0.14425f },
 124     { 0.25f, -0.14425f, -0.14425f,  0.14425f },
 125     { 0.25f,  0.14425f, -0.14425f, -0.14425f },
 126     { 0.25f, -0.14425f, -0.14425f, -0.14425f },
 127 };
 128 static ALfloat CubeEncoder[8][MAX_AMBI_COEFFS];
 129
 130 static alonce_flag encoder_inited = AL_ONCE_FLAG_INIT;
 131
 132 static void init_encoder(void)
 133 {
 134     CalcXYZCoeffs(-0.577350269f,  0.577350269f, -0.577350269f, CubeEncoder[0]);
 135     CalcXYZCoeffs( 0.577350269f,  0.577350269f, -0.577350269f, CubeEncoder[1]);
 136     CalcXYZCoeffs(-0.577350269f,  0.577350269f,  0.577350269f, CubeEncoder[2]);
 137     CalcXYZCoeffs( 0.577350269f,  0.577350269f,  0.577350269f, CubeEncoder[3]);
 138     CalcXYZCoeffs(-0.577350269f, -0.577350269f, -0.577350269f, CubeEncoder[4]);
 139     CalcXYZCoeffs( 0.577350269f, -0.577350269f, -0.577350269f, CubeEncoder[5]);
 140     CalcXYZCoeffs(-0.577350269f, -0.577350269f,  0.577350269f, CubeEncoder[6]);
 141     CalcXYZCoeffs( 0.577350269f, -0.577350269f,  0.577350269f, CubeEncoder[7]);
 142
 143     CalcXYZCoeffs(-0.707106781f,  0.0f, -0.707106781f, SquareEncoder[0]);
 144     CalcXYZCoeffs( 0.707106781f,  0.0f, -0.707106781f, SquareEncoder[1]);
 145     CalcXYZCoeffs(-0.707106781f,  0.0f,  0.707106781f, SquareEncoder[2]);
 146     CalcXYZCoeffs( 0.707106781f,  0.0f,  0.707106781f, SquareEncoder[3]);
 147 }
 148
 149
 150 /* NOTE: Low-frequency (LF) fields and BandSplitter filters are unused with
 151  * single-band decoding
 152  */
 153 typedef struct BFormatDec {
 154     alignas(16) ALfloat MatrixHF[MAX_OUTPUT_CHANNELS][MAX_AMBI_COEFFS];
 155     alignas(16) ALfloat MatrixLF[MAX_OUTPUT_CHANNELS][MAX_AMBI_COEFFS];
 156
 157     BandSplitter XOver[MAX_AMBI_COEFFS];
 158
 159     ALfloat (*Samples)[BUFFERSIZE];
 160     /* These two alias into Samples */
 161     ALfloat (*SamplesHF)[BUFFERSIZE];
 162     ALfloat (*SamplesLF)[BUFFERSIZE];
 163
 164     struct {
 165         const ALfloat (*restrict MatrixHF)[MAX_AMBI_COEFFS];
 166         const ALfloat (*restrict Encoder)[MAX_AMBI_COEFFS];
 167         ALuint NumChannels;
 168     } UpSampler;
 169
 170     ALuint NumChannels;
 171     ALboolean DualBand;
 172 } BFormatDec;
 173
 174 BFormatDec *bformatdec_alloc()
 175 {
 176     alcall_once(&encoder_inited, init_encoder);
 177     return al_calloc(16, sizeof(BFormatDec));
 178 }
 179
 180 void bformatdec_free(BFormatDec *dec)
 181 {
 182     if(dec)
 183     {
 184         al_free(dec->Samples);
 185         dec->Samples = NULL;
 186         dec->SamplesHF = NULL;
 187         dec->SamplesLF = NULL;
 188
 189         memset(dec, 0, sizeof(*dec));
 190         al_free(dec);
 191     }
 192 }
 193
 194 int bformatdec_getOrder(const struct BFormatDec *dec)
 195 {
 196     if(dec->NumChannels > 9) return 3;
 197     if(dec->NumChannels > 4) return 2;
 198     if(dec->NumChannels > 1) return 1;
 199     return 0;
 200 }
 201
 202 void bformatdec_reset(BFormatDec *dec, const AmbDecConf *conf, ALuint chancount, ALuint srate, const ALuint chanmap[MAX_OUTPUT_CHANNELS])
 203 {
 204     const ALfloat *coeff_scale = UnitScale;
 205     ALfloat ratio;
 206     ALuint i;
 207
 208     al_free(dec->Samples);
 209     dec->Samples = NULL;
 210     dec->SamplesHF = NULL;
 211     dec->SamplesLF = NULL;
 212
 213     dec->NumChannels = chancount;
 214     dec->Samples = al_calloc(16, dec->NumChannels * conf->FreqBands *
 215                                  sizeof(dec->Samples[0]));
 216     dec->SamplesHF = dec->Samples;
 217     dec->SamplesLF = dec->SamplesHF + dec->NumChannels;
 218
 219     if(conf->CoeffScale == ADS_SN3D)
 220         coeff_scale = SN3D2N3DScale;
 221     else if(conf->CoeffScale == ADS_FuMa)
 222         coeff_scale = FuMa2N3DScale;
 223
 224     if((conf->ChanMask & ~0x831b))
 225     {
 226         dec->UpSampler.MatrixHF = CubeMatrixHF;
 227         dec->UpSampler.Encoder = (const ALfloat(*)[MAX_AMBI_COEFFS])CubeEncoder;
 228         dec->UpSampler.NumChannels = 8;
 229     }
 230     else
 231     {
 232         dec->UpSampler.MatrixHF = SquareMatrixHF;
 233         dec->UpSampler.Encoder = (const ALfloat(*)[MAX_AMBI_COEFFS])SquareEncoder;
 234         dec->UpSampler.NumChannels = 4;
 235     }
 236
 237     if(conf->FreqBands == 1)
 238     {
 239         dec->DualBand = AL_FALSE;
 240         ratio = 1.0f;
 241     }
 242     else
 243     {
 244         dec->DualBand = AL_TRUE;
 245
 246         ratio = conf->XOverFreq / (ALfloat)srate;
 247         for(i = 0;i < MAX_AMBI_COEFFS;i++)
 248             bandsplit_init(&dec->XOver[i], ratio);
 249
 250         ratio = powf(10.0f, conf->XOverRatio / 40.0f);
 251         memset(dec->MatrixLF, 0, sizeof(dec->MatrixLF));
 252         for(i = 0;i < conf->NumSpeakers;i++)
 253         {
 254             ALuint chan = chanmap[i];
 255             ALuint j, k = 0;
 256             ALfloat gain;
 257
 258             for(j = 0;j < MAX_AMBI_COEFFS;j++)
 259             {
 260                 if(j == 0) gain = conf->LFOrderGain[0] / ratio;
 261                 else if(j == 1) gain = conf->LFOrderGain[1] / ratio;
 262                 else if(j == 4) gain = conf->LFOrderGain[2] / ratio;
 263                 else if(j == 9) gain = conf->LFOrderGain[3] / ratio;
 264                 if((conf->ChanMask&(1<<j)))
 265                     dec->MatrixLF[chan][j] = conf->LFMatrix[i][k++] / coeff_scale[j] * gain;
 266             }
 267         }
 268     }
 269
 270     memset(dec->MatrixHF, 0, sizeof(dec->MatrixHF));
 271     for(i = 0;i < conf->NumSpeakers;i++)
 272     {
 273         ALuint chan = chanmap[i];
 274         ALuint j, k = 0;
 275         ALfloat gain;
 276
 277         for(j = 0;j < MAX_AMBI_COEFFS;j++)
 278         {
 279             if(j == 0) gain = conf->HFOrderGain[0] * ratio;
 280             else if(j == 1) gain = conf->HFOrderGain[1] * ratio;
 281             else if(j == 4) gain = conf->HFOrderGain[2] * ratio;
 282             else if(j == 9) gain = conf->HFOrderGain[3] * ratio;
 283             if((conf->ChanMask&(1<<j)))
 284                 dec->MatrixHF[chan][j] = conf->HFMatrix[i][k++] / coeff_scale[j] * gain;
 285         }
 286     }
 287 }
 288
 289
 290 static void apply_row(ALfloat *out, const ALfloat *mtx, ALfloat (*restrict in)[BUFFERSIZE], ALuint inchans, ALuint todo)
 291 {
 292     ALuint c, i;
 293
 294     for(c = 0;c < inchans;c++)
 295     {
 296         ALfloat gain = mtx[c];
 297         if(!(fabsf(gain) > GAIN_SILENCE_THRESHOLD))
 298             continue;
 299         for(i = 0;i < todo;i++)
 300             out[i] += in[c][i] * gain;
 301     }
 302 }
 303
 304 void bformatdec_process(struct BFormatDec *dec, ALfloat (*restrict OutBuffer)[BUFFERSIZE], ALuint OutChannels, ALfloat (*restrict InSamples)[BUFFERSIZE], ALuint SamplesToDo)
 305 {
 306     ALuint chan, i;
 307
 308     if(dec->DualBand)
 309     {
 310         for(i = 0;i < dec->NumChannels;i++)
 311             bandsplit_process(&dec->XOver[i], dec->SamplesHF[i], dec->SamplesLF[i],
 312                               InSamples[i], SamplesToDo);
 313
 314         for(chan = 0;chan < OutChannels;chan++)
 315         {
 316             apply_row(OutBuffer[chan], dec->MatrixHF[chan], dec->SamplesHF,
 317                       dec->NumChannels, SamplesToDo);
 318             apply_row(OutBuffer[chan], dec->MatrixLF[chan], dec->SamplesLF,
 319                       dec->NumChannels, SamplesToDo);
 320         }
 321     }
 322     else
 323     {
 324         for(chan = 0;chan < OutChannels;chan++)
 325             apply_row(OutBuffer[chan], dec->MatrixHF[chan], InSamples,
 326                       dec->NumChannels, SamplesToDo);
 327     }
 328 }
 329
 330
 331 void bformatdec_upSample(struct BFormatDec *dec, ALfloat (*restrict OutBuffer)[BUFFERSIZE], ALfloat (*restrict InSamples)[BUFFERSIZE], ALuint InChannels, ALuint SamplesToDo)
 332 {
 333     ALuint i, j, k;
 334
 335     /* This up-sampler is very simplistic. It essentially decodes the first-
 336      * order content to a square channel array (or cube if height is desired),
 337      * then encodes those points onto the higher order soundfield.
 338      */
 339     for(k = 0;k < dec->UpSampler.NumChannels;k++)
 340     {
 341         memset(dec->Samples[0], 0, SamplesToDo*sizeof(ALfloat));
 342         apply_row(dec->Samples[0], dec->UpSampler.MatrixHF[k], InSamples,
 343                   InChannels, SamplesToDo);
 344
 345         for(j = 0;j < dec->NumChannels;j++)
 346         {
 347             ALfloat gain = dec->UpSampler.Encoder[k][j];
 348             if(!(fabsf(gain) > GAIN_SILENCE_THRESHOLD))
 349                 continue;
 350             for(i = 0;i < SamplesToDo;i++)
 351                 OutBuffer[j][i] += dec->Samples[0][i] * gain;
 352         }
 353     }
 354 }