nihav-codec-support/src/dsp/qmf.rs

   1 //! Quadrature Mirror Filter used in various audio codecs like MPEG Audio.
   2
   3 ///! 32-band QMF.
   4 pub struct QMF {
   5     hist:       [f32; 1024],
   6     pos:        usize,
   7 }
   8
   9 macro_rules! butterfly {
  10     (in; $src0: expr, $src1: expr, $dst0: expr, $dst1: expr, $scale: expr) => {
  11         $dst0 = $src0 + $src1;
  12         $dst1 = ($src0 - $src1) * $scale;
  13     };
  14     (rev; $a: expr, $b: expr) => {
  15         butterfly!(rev; $a, $b, std::f32::consts::FRAC_1_SQRT_2);
  16     };
  17     (rev; $a: expr, $b: expr, $scale: expr) => {
  18         let tmp = $a + $b;
  19         $b = ($a - $b) * $scale;
  20         $a = tmp;
  21     };
  22     (scal; $a: expr, $b: expr) => {
  23         butterfly!(scal; $a, $b, std::f32::consts::FRAC_1_SQRT_2);
  24     };
  25     (scal; $a: expr, $b: expr, $scale: expr) => {
  26         let tmp = $a + $b;
  27         $b = ($b - $a) * $scale;
  28         $a = tmp;
  29     };
  30 }
  31
  32 macro_rules! postadd {
  33     ($tmp: expr, $a0: expr, $b0: expr, $c0: expr, $d0: expr, $a1: expr, $b1: expr, $c1: expr, $d1: expr) => {
  34         $tmp[$c0] += $tmp[$d0];
  35         $tmp[$c1] += $tmp[$d1];
  36         $tmp[$a1] += $tmp[$c1];
  37         $tmp[$c1] += $tmp[$b1];
  38         $tmp[$b1] += $tmp[$d1];
  39     }
  40 }
  41
  42 fn dct32(src: &[f32; 32], dst: &mut [f32]) {
  43     let mut tmp = [0.0f32; 32];
  44
  45     butterfly!(in; src[ 0], src[31], tmp[ 0], tmp[31], 0.50060299823519627);
  46     butterfly!(in; src[ 1], src[30], tmp[ 1], tmp[30], 0.50547095989754365);
  47     butterfly!(in; src[ 2], src[29], tmp[ 2], tmp[29], 0.51544730992262455);
  48     butterfly!(in; src[ 3], src[28], tmp[ 3], tmp[28], 0.53104259108978413);
  49     butterfly!(in; src[ 4], src[27], tmp[ 4], tmp[27], 0.55310389603444454);
  50     butterfly!(in; src[ 5], src[26], tmp[ 5], tmp[26], 0.58293496820613389);
  51     butterfly!(in; src[ 6], src[25], tmp[ 6], tmp[25], 0.62250412303566482);
  52     butterfly!(in; src[ 7], src[24], tmp[ 7], tmp[24], 0.67480834145500568);
  53     butterfly!(in; src[ 8], src[23], tmp[ 8], tmp[23], 0.74453627100229858);
  54     butterfly!(in; src[ 9], src[22], tmp[ 9], tmp[22], 0.83934964541552681);
  55     butterfly!(in; src[10], src[21], tmp[10], tmp[21], 0.97256823786196078);
  56     butterfly!(in; src[11], src[20], tmp[11], tmp[20], 1.1694399334328847);
  57     butterfly!(in; src[12], src[19], tmp[12], tmp[19], 1.4841646163141662);
  58     butterfly!(in; src[13], src[18], tmp[13], tmp[18], 2.0577810099534108);
  59     butterfly!(in; src[14], src[17], tmp[14], tmp[17], 3.407608418468719);
  60     butterfly!(in; src[15], src[16], tmp[15], tmp[16], 10.190008123548033);
  61     butterfly!(rev;  tmp[ 0], tmp[15], 0.50241928618815568);
  62     butterfly!(rev;  tmp[ 1], tmp[14], 0.52249861493968885);
  63     butterfly!(rev;  tmp[ 2], tmp[13], 0.56694403481635769);
  64     butterfly!(rev;  tmp[ 3], tmp[12], 0.64682178335999008);
  65     butterfly!(rev;  tmp[ 4], tmp[11], 0.7881546234512502);
  66     butterfly!(rev;  tmp[ 5], tmp[10], 1.0606776859903471);
  67     butterfly!(rev;  tmp[ 6], tmp[ 9], 1.7224470982383342);
  68     butterfly!(rev;  tmp[ 7], tmp[ 8], 5.1011486186891553);
  69     butterfly!(scal; tmp[16], tmp[31], 0.50241928618815568);
  70     butterfly!(scal; tmp[17], tmp[30], 0.52249861493968885);
  71     butterfly!(scal; tmp[18], tmp[29], 0.56694403481635769);
  72     butterfly!(scal; tmp[19], tmp[28], 0.64682178335999008);
  73     butterfly!(scal; tmp[20], tmp[27], 0.7881546234512502);
  74     butterfly!(scal; tmp[21], tmp[26], 1.0606776859903471);
  75     butterfly!(scal; tmp[22], tmp[25], 1.7224470982383342);
  76     butterfly!(scal; tmp[23], tmp[24], 5.1011486186891553);
  77     butterfly!(rev;  tmp[ 0], tmp[ 7], 0.50979557910415918);
  78     butterfly!(rev;  tmp[ 1], tmp[ 6], 0.60134488693504529);
  79     butterfly!(rev;  tmp[ 2], tmp[ 5], 0.89997622313641557);
  80     butterfly!(rev;  tmp[ 3], tmp[ 4], 2.5629154477415055);
  81     butterfly!(rev;  tmp[16], tmp[23], 0.50979557910415918);
  82     butterfly!(rev;  tmp[17], tmp[22], 0.60134488693504529);
  83     butterfly!(rev;  tmp[18], tmp[21], 0.89997622313641557);
  84     butterfly!(rev;  tmp[19], tmp[20], 2.5629154477415055);
  85     butterfly!(scal; tmp[ 8], tmp[15], 0.50979557910415918);
  86     butterfly!(scal; tmp[ 9], tmp[14], 0.60134488693504529);
  87     butterfly!(scal; tmp[10], tmp[13], 0.89997622313641557);
  88     butterfly!(scal; tmp[11], tmp[12], 2.5629154477415055);
  89     butterfly!(scal; tmp[24], tmp[31], 0.50979557910415918);
  90     butterfly!(scal; tmp[25], tmp[30], 0.60134488693504529);
  91     butterfly!(scal; tmp[26], tmp[29], 0.89997622313641557);
  92     butterfly!(scal; tmp[27], tmp[28], 2.5629154477415055);
  93     butterfly!(rev;  tmp[ 0], tmp[ 3], 0.54119610014619701);
  94     butterfly!(rev;  tmp[ 1], tmp[ 2], 1.3065629648763764);
  95     butterfly!(rev;  tmp[ 8], tmp[11], 0.54119610014619701);
  96     butterfly!(rev;  tmp[ 9], tmp[10], 1.3065629648763764);
  97     butterfly!(rev;  tmp[16], tmp[19], 0.54119610014619701);
  98     butterfly!(rev;  tmp[17], tmp[18], 1.3065629648763764);
  99     butterfly!(rev;  tmp[24], tmp[27], 0.54119610014619701);
 100     butterfly!(rev;  tmp[25], tmp[26], 1.3065629648763764);
 101     butterfly!(scal; tmp[ 4], tmp[ 7], 0.54119610014619701);
 102     butterfly!(scal; tmp[ 5], tmp[ 6], 1.3065629648763764);
 103     butterfly!(scal; tmp[12], tmp[15], 0.54119610014619701);
 104     butterfly!(scal; tmp[13], tmp[14], 1.3065629648763764);
 105     butterfly!(scal; tmp[20], tmp[23], 0.54119610014619701);
 106     butterfly!(scal; tmp[21], tmp[22], 1.3065629648763764);
 107     butterfly!(scal; tmp[28], tmp[31], 0.54119610014619701);
 108     butterfly!(scal; tmp[29], tmp[30], 1.3065629648763764);
 109     butterfly!(rev;  tmp[ 0], tmp[ 1]);
 110     butterfly!(rev;  tmp[ 4], tmp[ 5]);
 111     butterfly!(rev;  tmp[ 8], tmp[ 9]);
 112     butterfly!(rev;  tmp[12], tmp[13]);
 113     butterfly!(rev;  tmp[16], tmp[17]);
 114     butterfly!(rev;  tmp[20], tmp[21]);
 115     butterfly!(rev;  tmp[24], tmp[25]);
 116     butterfly!(rev;  tmp[28], tmp[29]);
 117     butterfly!(scal; tmp[ 2], tmp[ 3]);
 118     butterfly!(scal; tmp[ 6], tmp[ 7]);
 119     butterfly!(scal; tmp[10], tmp[11]);
 120     butterfly!(scal; tmp[14], tmp[15]);
 121     butterfly!(scal; tmp[18], tmp[19]);
 122     butterfly!(scal; tmp[22], tmp[23]);
 123     butterfly!(scal; tmp[26], tmp[27]);
 124     butterfly!(scal; tmp[30], tmp[31]);
 125
 126     postadd!(tmp, 0,  1,  2,  3,  4,  5,  6,  7);
 127     postadd!(tmp, 8,  9, 10, 11, 12, 13, 14, 15);
 128     postadd!(tmp,16, 17, 18, 19, 20, 21, 22, 23);
 129     postadd!(tmp,24, 25, 26, 27, 28, 29, 30, 31);
 130
 131     dst[ 0] = tmp[0];
 132     dst[16] = tmp[1];
 133     dst[ 8] = tmp[2];
 134     dst[24] = tmp[3];
 135     dst[ 4] = tmp[4];
 136     dst[20] = tmp[5];
 137     dst[12] = tmp[6];
 138     dst[28] = tmp[7];
 139
 140     dst[ 2] = tmp[ 8] + tmp[12];
 141     dst[18] = tmp[ 9] + tmp[13];
 142     dst[10] = tmp[10] + tmp[14];
 143     dst[26] = tmp[11] + tmp[15];
 144
 145     dst[ 6] = tmp[12] + tmp[10];
 146     dst[22] = tmp[13] + tmp[11];
 147
 148     dst[14] = tmp[14] + tmp[ 9];
 149
 150     dst[30] = tmp[15];
 151
 152     tmp[24] += tmp[28];
 153     tmp[28] += tmp[26];
 154     tmp[26] += tmp[30];
 155     tmp[30] += tmp[25];
 156     tmp[25] += tmp[29];
 157     tmp[29] += tmp[27];
 158     tmp[27] += tmp[31];
 159
 160     dst[ 1] = tmp[16] + tmp[24];
 161     dst[17] = tmp[17] + tmp[25];
 162     dst[ 9] = tmp[18] + tmp[26];
 163     dst[25] = tmp[19] + tmp[27];
 164     dst[ 5] = tmp[20] + tmp[28];
 165     dst[21] = tmp[21] + tmp[29];
 166     dst[13] = tmp[22] + tmp[30];
 167     dst[29] = tmp[23] + tmp[31];
 168
 169     dst[ 3] = tmp[24] + tmp[20];
 170     dst[19] = tmp[25] + tmp[21];
 171     dst[11] = tmp[26] + tmp[22];
 172     dst[27] = tmp[27] + tmp[23];
 173
 174     dst[ 7] = tmp[28] + tmp[18];
 175     dst[23] = tmp[29] + tmp[19];
 176
 177     dst[15] = tmp[30] + tmp[17];
 178
 179     dst[31] = tmp[31];
 180 }
 181
 182 impl QMF {
 183     ///! Constructs a new `QMF` context.
 184     pub fn new() -> Self {
 185         Self {
 186             hist:   [0.0; 1024],
 187             pos:    0,
 188         }
 189     }
 190     ///! Performs reconstruction of audio from provided sub-band samples.
 191     pub fn synth(&mut self, src: &[f32; 32], dst: &mut [f32]) {
 192         self.pos = self.pos.wrapping_sub(32) & 0x1FF;
 193         dct32(src, &mut self.hist[self.pos..][..32]);
 194
 195         for i in 0..16 {
 196             let mut acc = 0.0;
 197             for j in (0..512).step_by(64) {
 198                 acc += QMF_WINDOW[i + j]      * self.hist[(j + 16      + i + self.pos) & 0x1FF];
 199                 acc -= QMF_WINDOW[i + j + 32] * self.hist[(j + 16 + 32 - i + self.pos) & 0x1FF];
 200             }
 201             dst[i] = acc;
 202         }
 203         let mut acc = 0.0;
 204         for j in (0..512).step_by(64) {
 205             acc -= QMF_WINDOW[j + 32 + 16] * self.hist[(j + 32 + self.pos) & 0x1FF];
 206         }
 207         dst[16] = acc;
 208         for i in 17..32 {
 209             let mut acc = 0.0;
 210             for j in (0..512).step_by(64) {
 211                 acc -= QMF_WINDOW[i + j]      * self.hist[(j + 32 + 16 - i + self.pos) & 0x1FF];
 212                 acc -= QMF_WINDOW[i + j + 32] * self.hist[(j + 32 - 16 + i + self.pos) & 0x1FF];
 213             }
 214             dst[i] = acc;
 215         }
 216     }
 217     ///! Resets internal filter state.
 218     pub fn reset(&mut self) {
 219         self.hist = [0.0; 1024];
 220         self.pos  = 0;
 221     }
 222 }
 223
 224 impl Default for QMF {
 225     fn default() -> Self {
 226         Self::new()
 227     }
 228 }
 229
 230 const QMF_WINDOW: [f32; 512] = [
 231   0.000000000, -0.000015259, -0.000015259, -0.000015259,
 232  -0.000015259, -0.000015259, -0.000015259, -0.000030518,
 233  -0.000030518, -0.000030518, -0.000030518, -0.000045776,
 234  -0.000045776, -0.000061035, -0.000061035, -0.000076294,
 235  -0.000076294, -0.000091553, -0.000106812, -0.000106812,
 236  -0.000122070, -0.000137329, -0.000152588, -0.000167847,
 237  -0.000198364, -0.000213623, -0.000244141, -0.000259399,
 238  -0.000289917, -0.000320435, -0.000366211, -0.000396729,
 239  -0.000442505, -0.000473022, -0.000534058, -0.000579834,
 240  -0.000625610, -0.000686646, -0.000747681, -0.000808716,
 241  -0.000885010, -0.000961304, -0.001037598, -0.001113892,
 242  -0.001205444, -0.001296997, -0.001388550, -0.001480103,
 243  -0.001586914, -0.001693726, -0.001785278, -0.001907349,
 244  -0.002014160, -0.002120972, -0.002243042, -0.002349854,
 245  -0.002456665, -0.002578735, -0.002685547, -0.002792358,
 246  -0.002899170, -0.002990723, -0.003082275, -0.003173828,
 247   0.003250122,  0.003326416,  0.003387451,  0.003433228,
 248   0.003463745,  0.003479004,  0.003479004,  0.003463745,
 249   0.003417969,  0.003372192,  0.003280640,  0.003173828,
 250   0.003051758,  0.002883911,  0.002700806,  0.002487183,
 251   0.002227783,  0.001937866,  0.001617432,  0.001266479,
 252   0.000869751,  0.000442505, -0.000030518, -0.000549316,
 253  -0.001098633, -0.001693726, -0.002334595, -0.003005981,
 254  -0.003723145, -0.004486084, -0.005294800, -0.006118774,
 255  -0.007003784, -0.007919312, -0.008865356, -0.009841919,
 256  -0.010848999, -0.011886597, -0.012939453, -0.014022827,
 257  -0.015121460, -0.016235352, -0.017349243, -0.018463135,
 258  -0.019577026, -0.020690918, -0.021789551, -0.022857666,
 259  -0.023910522, -0.024932861, -0.025909424, -0.026840210,
 260  -0.027725220, -0.028533936, -0.029281616, -0.029937744,
 261  -0.030532837, -0.031005859, -0.031387329, -0.031661987,
 262  -0.031814575, -0.031845093, -0.031738281, -0.031478882,
 263   0.031082153,  0.030517578,  0.029785156,  0.028884888,
 264   0.027801514,  0.026535034,  0.025085449,  0.023422241,
 265   0.021575928,  0.019531250,  0.017257690,  0.014801025,
 266   0.012115479,  0.009231567,  0.006134033,  0.002822876,
 267  -0.000686646, -0.004394531, -0.008316040, -0.012420654,
 268  -0.016708374, -0.021179199, -0.025817871, -0.030609131,
 269  -0.035552979, -0.040634155, -0.045837402, -0.051132202,
 270  -0.056533813, -0.061996460, -0.067520142, -0.073059082,
 271  -0.078628540, -0.084182739, -0.089706421, -0.095169067,
 272  -0.100540161, -0.105819702, -0.110946655, -0.115921021,
 273  -0.120697021, -0.125259399, -0.129562378, -0.133590698,
 274  -0.137298584, -0.140670776, -0.143676758, -0.146255493,
 275  -0.148422241, -0.150115967, -0.151306152, -0.151962280,
 276  -0.152069092, -0.151596069, -0.150497437, -0.148773193,
 277  -0.146362305, -0.143264771, -0.139450073, -0.134887695,
 278  -0.129577637, -0.123474121, -0.116577148, -0.108856201,
 279   0.100311279,  0.090927124,  0.080688477,  0.069595337,
 280   0.057617187,  0.044784546,  0.031082153,  0.016510010,
 281   0.001068115, -0.015228271, -0.032379150, -0.050354004,
 282  -0.069168091, -0.088775635, -0.109161377, -0.130310059,
 283  -0.152206421, -0.174789429, -0.198059082, -0.221984863,
 284  -0.246505737, -0.271591187, -0.297210693, -0.323318481,
 285  -0.349868774, -0.376800537, -0.404083252, -0.431655884,
 286  -0.459472656, -0.487472534, -0.515609741, -0.543823242,
 287  -0.572036743, -0.600219727, -0.628295898, -0.656219482,
 288  -0.683914185, -0.711318970, -0.738372803, -0.765029907,
 289  -0.791213989, -0.816864014, -0.841949463, -0.866363525,
 290  -0.890090942, -0.913055420, -0.935195923, -0.956481934,
 291  -0.976852417, -0.996246338, -1.014617920, -1.031936646,
 292  -1.048156738, -1.063217163, -1.077117920, -1.089782715,
 293  -1.101211548, -1.111373901, -1.120223999, -1.127746582,
 294  -1.133926392, -1.138763428, -1.142211914, -1.144287109,
 295   1.144989014,  1.144287109,  1.142211914,  1.138763428,
 296   1.133926392,  1.127746582,  1.120223999,  1.111373901,
 297   1.101211548,  1.089782715,  1.077117920,  1.063217163,
 298   1.048156738,  1.031936646,  1.014617920,  0.996246338,
 299   0.976852417,  0.956481934,  0.935195923,  0.913055420,
 300   0.890090942,  0.866363525,  0.841949463,  0.816864014,
 301   0.791213989,  0.765029907,  0.738372803,  0.711318970,
 302   0.683914185,  0.656219482,  0.628295898,  0.600219727,
 303   0.572036743,  0.543823242,  0.515609741,  0.487472534,
 304   0.459472656,  0.431655884,  0.404083252,  0.376800537,
 305   0.349868774,  0.323318481,  0.297210693,  0.271591187,
 306   0.246505737,  0.221984863,  0.198059082,  0.174789429,
 307   0.152206421,  0.130310059,  0.109161377,  0.088775635,
 308   0.069168091,  0.050354004,  0.032379150,  0.015228271,
 309  -0.001068115, -0.016510010, -0.031082153, -0.044784546,
 310  -0.057617187, -0.069595337, -0.080688477, -0.090927124,
 311   0.100311279,  0.108856201,  0.116577148,  0.123474121,
 312   0.129577637,  0.134887695,  0.139450073,  0.143264771,
 313   0.146362305,  0.148773193,  0.150497437,  0.151596069,
 314   0.152069092,  0.151962280,  0.151306152,  0.150115967,
 315   0.148422241,  0.146255493,  0.143676758,  0.140670776,
 316   0.137298584,  0.133590698,  0.129562378,  0.125259399,
 317   0.120697021,  0.115921021,  0.110946655,  0.105819702,
 318   0.100540161,  0.095169067,  0.089706421,  0.084182739,
 319   0.078628540,  0.073059082,  0.067520142,  0.061996460,
 320   0.056533813,  0.051132202,  0.045837402,  0.040634155,
 321   0.035552979,  0.030609131,  0.025817871,  0.021179199,
 322   0.016708374,  0.012420654,  0.008316040,  0.004394531,
 323   0.000686646, -0.002822876, -0.006134033, -0.009231567,
 324  -0.012115479, -0.014801025, -0.017257690, -0.019531250,
 325  -0.021575928, -0.023422241, -0.025085449, -0.026535034,
 326  -0.027801514, -0.028884888, -0.029785156, -0.030517578,
 327   0.031082153,  0.031478882,  0.031738281,  0.031845093,
 328   0.031814575,  0.031661987,  0.031387329,  0.031005859,
 329   0.030532837,  0.029937744,  0.029281616,  0.028533936,
 330   0.027725220,  0.026840210,  0.025909424,  0.024932861,
 331   0.023910522,  0.022857666,  0.021789551,  0.020690918,
 332   0.019577026,  0.018463135,  0.017349243,  0.016235352,
 333   0.015121460,  0.014022827,  0.012939453,  0.011886597,
 334   0.010848999,  0.009841919,  0.008865356,  0.007919312,
 335   0.007003784,  0.006118774,  0.005294800,  0.004486084,
 336   0.003723145,  0.003005981,  0.002334595,  0.001693726,
 337   0.001098633,  0.000549316,  0.000030518, -0.000442505,
 338  -0.000869751, -0.001266479, -0.001617432, -0.001937866,
 339  -0.002227783, -0.002487183, -0.002700806, -0.002883911,
 340  -0.003051758, -0.003173828, -0.003280640, -0.003372192,
 341  -0.003417969, -0.003463745, -0.003479004, -0.003479004,
 342  -0.003463745, -0.003433228, -0.003387451, -0.003326416,
 343   0.003250122,  0.003173828,  0.003082275,  0.002990723,
 344   0.002899170,  0.002792358,  0.002685547,  0.002578735,
 345   0.002456665,  0.002349854,  0.002243042,  0.002120972,
 346   0.002014160,  0.001907349,  0.001785278,  0.001693726,
 347   0.001586914,  0.001480103,  0.001388550,  0.001296997,
 348   0.001205444,  0.001113892,  0.001037598,  0.000961304,
 349   0.000885010,  0.000808716,  0.000747681,  0.000686646,
 350   0.000625610,  0.000579834,  0.000534058,  0.000473022,
 351   0.000442505,  0.000396729,  0.000366211,  0.000320435,
 352   0.000289917,  0.000259399,  0.000244141,  0.000213623,
 353   0.000198364,  0.000167847,  0.000152588,  0.000137329,
 354   0.000122070,  0.000106812,  0.000106812,  0.000091553,
 355   0.000076294,  0.000076294,  0.000061035,  0.000061035,
 356   0.000045776,  0.000045776,  0.000030518,  0.000030518,
 357   0.000030518,  0.000030518,  0.000015259,  0.000015259,
 358   0.000015259,  0.000015259,  0.000015259,  0.000015259,
 359 ];