Alexandre Lision | 8538238 | 2014-01-27 15:54:16 -0500 | [diff] [blame] | 1 | /* Copyright (c) 2011 Xiph.Org Foundation |
| 2 | Written by Jean-Marc Valin */ |
| 3 | /* |
| 4 | Redistribution and use in source and binary forms, with or without |
| 5 | modification, are permitted provided that the following conditions |
| 6 | are met: |
| 7 | |
| 8 | - Redistributions of source code must retain the above copyright |
| 9 | notice, this list of conditions and the following disclaimer. |
| 10 | |
| 11 | - Redistributions in binary form must reproduce the above copyright |
| 12 | notice, this list of conditions and the following disclaimer in the |
| 13 | documentation and/or other materials provided with the distribution. |
| 14 | |
| 15 | THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS |
| 16 | ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT |
| 17 | LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR |
| 18 | A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE FOUNDATION OR |
| 19 | CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, |
| 20 | EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, |
| 21 | PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR |
| 22 | PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF |
| 23 | LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING |
| 24 | NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS |
| 25 | SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. |
| 26 | */ |
| 27 | |
| 28 | #ifdef HAVE_CONFIG_H |
| 29 | #include "config.h" |
| 30 | #endif |
| 31 | |
| 32 | #include "kiss_fft.h" |
| 33 | #include "celt.h" |
| 34 | #include "modes.h" |
| 35 | #include "arch.h" |
| 36 | #include "quant_bands.h" |
| 37 | #include <stdio.h> |
| 38 | #include "analysis.h" |
| 39 | #include "mlp.h" |
| 40 | #include "stack_alloc.h" |
| 41 | |
| 42 | extern const MLP net; |
| 43 | |
| 44 | #ifndef M_PI |
| 45 | #define M_PI 3.141592653 |
| 46 | #endif |
| 47 | |
| 48 | static const float dct_table[128] = { |
| 49 | 0.250000f, 0.250000f, 0.250000f, 0.250000f, 0.250000f, 0.250000f, 0.250000f, 0.250000f, |
| 50 | 0.250000f, 0.250000f, 0.250000f, 0.250000f, 0.250000f, 0.250000f, 0.250000f, 0.250000f, |
| 51 | 0.351851f, 0.338330f, 0.311806f, 0.273300f, 0.224292f, 0.166664f, 0.102631f, 0.034654f, |
| 52 | -0.034654f,-0.102631f,-0.166664f,-0.224292f,-0.273300f,-0.311806f,-0.338330f,-0.351851f, |
| 53 | 0.346760f, 0.293969f, 0.196424f, 0.068975f,-0.068975f,-0.196424f,-0.293969f,-0.346760f, |
| 54 | -0.346760f,-0.293969f,-0.196424f,-0.068975f, 0.068975f, 0.196424f, 0.293969f, 0.346760f, |
| 55 | 0.338330f, 0.224292f, 0.034654f,-0.166664f,-0.311806f,-0.351851f,-0.273300f,-0.102631f, |
| 56 | 0.102631f, 0.273300f, 0.351851f, 0.311806f, 0.166664f,-0.034654f,-0.224292f,-0.338330f, |
| 57 | 0.326641f, 0.135299f,-0.135299f,-0.326641f,-0.326641f,-0.135299f, 0.135299f, 0.326641f, |
| 58 | 0.326641f, 0.135299f,-0.135299f,-0.326641f,-0.326641f,-0.135299f, 0.135299f, 0.326641f, |
| 59 | 0.311806f, 0.034654f,-0.273300f,-0.338330f,-0.102631f, 0.224292f, 0.351851f, 0.166664f, |
| 60 | -0.166664f,-0.351851f,-0.224292f, 0.102631f, 0.338330f, 0.273300f,-0.034654f,-0.311806f, |
| 61 | 0.293969f,-0.068975f,-0.346760f,-0.196424f, 0.196424f, 0.346760f, 0.068975f,-0.293969f, |
| 62 | -0.293969f, 0.068975f, 0.346760f, 0.196424f,-0.196424f,-0.346760f,-0.068975f, 0.293969f, |
| 63 | 0.273300f,-0.166664f,-0.338330f, 0.034654f, 0.351851f, 0.102631f,-0.311806f,-0.224292f, |
| 64 | 0.224292f, 0.311806f,-0.102631f,-0.351851f,-0.034654f, 0.338330f, 0.166664f,-0.273300f, |
| 65 | }; |
| 66 | |
| 67 | static const float analysis_window[240] = { |
| 68 | 0.000043f, 0.000171f, 0.000385f, 0.000685f, 0.001071f, 0.001541f, 0.002098f, 0.002739f, |
| 69 | 0.003466f, 0.004278f, 0.005174f, 0.006156f, 0.007222f, 0.008373f, 0.009607f, 0.010926f, |
| 70 | 0.012329f, 0.013815f, 0.015385f, 0.017037f, 0.018772f, 0.020590f, 0.022490f, 0.024472f, |
| 71 | 0.026535f, 0.028679f, 0.030904f, 0.033210f, 0.035595f, 0.038060f, 0.040604f, 0.043227f, |
| 72 | 0.045928f, 0.048707f, 0.051564f, 0.054497f, 0.057506f, 0.060591f, 0.063752f, 0.066987f, |
| 73 | 0.070297f, 0.073680f, 0.077136f, 0.080665f, 0.084265f, 0.087937f, 0.091679f, 0.095492f, |
| 74 | 0.099373f, 0.103323f, 0.107342f, 0.111427f, 0.115579f, 0.119797f, 0.124080f, 0.128428f, |
| 75 | 0.132839f, 0.137313f, 0.141849f, 0.146447f, 0.151105f, 0.155823f, 0.160600f, 0.165435f, |
| 76 | 0.170327f, 0.175276f, 0.180280f, 0.185340f, 0.190453f, 0.195619f, 0.200838f, 0.206107f, |
| 77 | 0.211427f, 0.216797f, 0.222215f, 0.227680f, 0.233193f, 0.238751f, 0.244353f, 0.250000f, |
| 78 | 0.255689f, 0.261421f, 0.267193f, 0.273005f, 0.278856f, 0.284744f, 0.290670f, 0.296632f, |
| 79 | 0.302628f, 0.308658f, 0.314721f, 0.320816f, 0.326941f, 0.333097f, 0.339280f, 0.345492f, |
| 80 | 0.351729f, 0.357992f, 0.364280f, 0.370590f, 0.376923f, 0.383277f, 0.389651f, 0.396044f, |
| 81 | 0.402455f, 0.408882f, 0.415325f, 0.421783f, 0.428254f, 0.434737f, 0.441231f, 0.447736f, |
| 82 | 0.454249f, 0.460770f, 0.467298f, 0.473832f, 0.480370f, 0.486912f, 0.493455f, 0.500000f, |
| 83 | 0.506545f, 0.513088f, 0.519630f, 0.526168f, 0.532702f, 0.539230f, 0.545751f, 0.552264f, |
| 84 | 0.558769f, 0.565263f, 0.571746f, 0.578217f, 0.584675f, 0.591118f, 0.597545f, 0.603956f, |
| 85 | 0.610349f, 0.616723f, 0.623077f, 0.629410f, 0.635720f, 0.642008f, 0.648271f, 0.654508f, |
| 86 | 0.660720f, 0.666903f, 0.673059f, 0.679184f, 0.685279f, 0.691342f, 0.697372f, 0.703368f, |
| 87 | 0.709330f, 0.715256f, 0.721144f, 0.726995f, 0.732807f, 0.738579f, 0.744311f, 0.750000f, |
| 88 | 0.755647f, 0.761249f, 0.766807f, 0.772320f, 0.777785f, 0.783203f, 0.788573f, 0.793893f, |
| 89 | 0.799162f, 0.804381f, 0.809547f, 0.814660f, 0.819720f, 0.824724f, 0.829673f, 0.834565f, |
| 90 | 0.839400f, 0.844177f, 0.848895f, 0.853553f, 0.858151f, 0.862687f, 0.867161f, 0.871572f, |
| 91 | 0.875920f, 0.880203f, 0.884421f, 0.888573f, 0.892658f, 0.896677f, 0.900627f, 0.904508f, |
| 92 | 0.908321f, 0.912063f, 0.915735f, 0.919335f, 0.922864f, 0.926320f, 0.929703f, 0.933013f, |
| 93 | 0.936248f, 0.939409f, 0.942494f, 0.945503f, 0.948436f, 0.951293f, 0.954072f, 0.956773f, |
| 94 | 0.959396f, 0.961940f, 0.964405f, 0.966790f, 0.969096f, 0.971321f, 0.973465f, 0.975528f, |
| 95 | 0.977510f, 0.979410f, 0.981228f, 0.982963f, 0.984615f, 0.986185f, 0.987671f, 0.989074f, |
| 96 | 0.990393f, 0.991627f, 0.992778f, 0.993844f, 0.994826f, 0.995722f, 0.996534f, 0.997261f, |
| 97 | 0.997902f, 0.998459f, 0.998929f, 0.999315f, 0.999615f, 0.999829f, 0.999957f, 1.000000f, |
| 98 | }; |
| 99 | |
| 100 | static const int tbands[NB_TBANDS+1] = { |
| 101 | 2, 4, 6, 8, 10, 12, 14, 16, 20, 24, 28, 32, 40, 48, 56, 68, 80, 96, 120 |
| 102 | }; |
| 103 | |
| 104 | static const int extra_bands[NB_TOT_BANDS+1] = { |
| 105 | 1, 2, 4, 6, 8, 10, 12, 14, 16, 20, 24, 28, 32, 40, 48, 56, 68, 80, 96, 120, 160, 200 |
| 106 | }; |
| 107 | |
| 108 | /*static const float tweight[NB_TBANDS+1] = { |
| 109 | .3, .4, .5, .6, .7, .8, .9, 1., 1., 1., 1., 1., 1., 1., .8, .7, .6, .5 |
| 110 | };*/ |
| 111 | |
| 112 | #define NB_TONAL_SKIP_BANDS 9 |
| 113 | |
| 114 | #define cA 0.43157974f |
| 115 | #define cB 0.67848403f |
| 116 | #define cC 0.08595542f |
| 117 | #define cE ((float)M_PI/2) |
| 118 | static OPUS_INLINE float fast_atan2f(float y, float x) { |
| 119 | float x2, y2; |
| 120 | /* Should avoid underflow on the values we'll get */ |
| 121 | if (ABS16(x)+ABS16(y)<1e-9f) |
| 122 | { |
| 123 | x*=1e12f; |
| 124 | y*=1e12f; |
| 125 | } |
| 126 | x2 = x*x; |
| 127 | y2 = y*y; |
| 128 | if(x2<y2){ |
| 129 | float den = (y2 + cB*x2) * (y2 + cC*x2); |
| 130 | if (den!=0) |
| 131 | return -x*y*(y2 + cA*x2) / den + (y<0 ? -cE : cE); |
| 132 | else |
| 133 | return (y<0 ? -cE : cE); |
| 134 | }else{ |
| 135 | float den = (x2 + cB*y2) * (x2 + cC*y2); |
| 136 | if (den!=0) |
| 137 | return x*y*(x2 + cA*y2) / den + (y<0 ? -cE : cE) - (x*y<0 ? -cE : cE); |
| 138 | else |
| 139 | return (y<0 ? -cE : cE) - (x*y<0 ? -cE : cE); |
| 140 | } |
| 141 | } |
| 142 | |
| 143 | void tonality_get_info(TonalityAnalysisState *tonal, AnalysisInfo *info_out, int len) |
| 144 | { |
| 145 | int pos; |
| 146 | int curr_lookahead; |
| 147 | float psum; |
| 148 | int i; |
| 149 | |
| 150 | pos = tonal->read_pos; |
| 151 | curr_lookahead = tonal->write_pos-tonal->read_pos; |
| 152 | if (curr_lookahead<0) |
| 153 | curr_lookahead += DETECT_SIZE; |
| 154 | |
| 155 | if (len > 480 && pos != tonal->write_pos) |
| 156 | { |
| 157 | pos++; |
| 158 | if (pos==DETECT_SIZE) |
| 159 | pos=0; |
| 160 | } |
| 161 | if (pos == tonal->write_pos) |
| 162 | pos--; |
| 163 | if (pos<0) |
| 164 | pos = DETECT_SIZE-1; |
| 165 | OPUS_COPY(info_out, &tonal->info[pos], 1); |
| 166 | tonal->read_subframe += len/120; |
| 167 | while (tonal->read_subframe>=4) |
| 168 | { |
| 169 | tonal->read_subframe -= 4; |
| 170 | tonal->read_pos++; |
| 171 | } |
| 172 | if (tonal->read_pos>=DETECT_SIZE) |
| 173 | tonal->read_pos-=DETECT_SIZE; |
| 174 | |
| 175 | /* Compensate for the delay in the features themselves. |
| 176 | FIXME: Need a better estimate the 10 I just made up */ |
| 177 | curr_lookahead = IMAX(curr_lookahead-10, 0); |
| 178 | |
| 179 | psum=0; |
| 180 | /* Summing the probability of transition patterns that involve music at |
| 181 | time (DETECT_SIZE-curr_lookahead-1) */ |
| 182 | for (i=0;i<DETECT_SIZE-curr_lookahead;i++) |
| 183 | psum += tonal->pmusic[i]; |
| 184 | for (;i<DETECT_SIZE;i++) |
| 185 | psum += tonal->pspeech[i]; |
| 186 | psum = psum*tonal->music_confidence + (1-psum)*tonal->speech_confidence; |
| 187 | /*printf("%f %f %f\n", psum, info_out->music_prob, info_out->tonality);*/ |
| 188 | |
| 189 | info_out->music_prob = psum; |
| 190 | } |
| 191 | |
| 192 | void tonality_analysis(TonalityAnalysisState *tonal, AnalysisInfo *info_out, const CELTMode *celt_mode, const void *x, int len, int offset, int c1, int c2, int C, int lsb_depth, downmix_func downmix) |
| 193 | { |
| 194 | int i, b; |
| 195 | const kiss_fft_state *kfft; |
| 196 | VARDECL(kiss_fft_cpx, in); |
| 197 | VARDECL(kiss_fft_cpx, out); |
| 198 | int N = 480, N2=240; |
| 199 | float * OPUS_RESTRICT A = tonal->angle; |
| 200 | float * OPUS_RESTRICT dA = tonal->d_angle; |
| 201 | float * OPUS_RESTRICT d2A = tonal->d2_angle; |
| 202 | VARDECL(float, tonality); |
| 203 | VARDECL(float, noisiness); |
| 204 | float band_tonality[NB_TBANDS]; |
| 205 | float logE[NB_TBANDS]; |
| 206 | float BFCC[8]; |
| 207 | float features[25]; |
| 208 | float frame_tonality; |
| 209 | float max_frame_tonality; |
| 210 | /*float tw_sum=0;*/ |
| 211 | float frame_noisiness; |
| 212 | const float pi4 = (float)(M_PI*M_PI*M_PI*M_PI); |
| 213 | float slope=0; |
| 214 | float frame_stationarity; |
| 215 | float relativeE; |
| 216 | float frame_probs[2]; |
| 217 | float alpha, alphaE, alphaE2; |
| 218 | float frame_loudness; |
| 219 | float bandwidth_mask; |
| 220 | int bandwidth=0; |
| 221 | float maxE = 0; |
| 222 | float noise_floor; |
| 223 | int remaining; |
| 224 | AnalysisInfo *info; |
| 225 | SAVE_STACK; |
| 226 | |
| 227 | tonal->last_transition++; |
| 228 | alpha = 1.f/IMIN(20, 1+tonal->count); |
| 229 | alphaE = 1.f/IMIN(50, 1+tonal->count); |
| 230 | alphaE2 = 1.f/IMIN(1000, 1+tonal->count); |
| 231 | |
| 232 | if (tonal->count<4) |
| 233 | tonal->music_prob = .5; |
| 234 | kfft = celt_mode->mdct.kfft[0]; |
| 235 | if (tonal->count==0) |
| 236 | tonal->mem_fill = 240; |
| 237 | downmix(x, &tonal->inmem[tonal->mem_fill], IMIN(len, ANALYSIS_BUF_SIZE-tonal->mem_fill), offset, c1, c2, C); |
| 238 | if (tonal->mem_fill+len < ANALYSIS_BUF_SIZE) |
| 239 | { |
| 240 | tonal->mem_fill += len; |
| 241 | /* Don't have enough to update the analysis */ |
| 242 | RESTORE_STACK; |
| 243 | return; |
| 244 | } |
| 245 | info = &tonal->info[tonal->write_pos++]; |
| 246 | if (tonal->write_pos>=DETECT_SIZE) |
| 247 | tonal->write_pos-=DETECT_SIZE; |
| 248 | |
| 249 | ALLOC(in, 480, kiss_fft_cpx); |
| 250 | ALLOC(out, 480, kiss_fft_cpx); |
| 251 | ALLOC(tonality, 240, float); |
| 252 | ALLOC(noisiness, 240, float); |
| 253 | for (i=0;i<N2;i++) |
| 254 | { |
| 255 | float w = analysis_window[i]; |
| 256 | in[i].r = (kiss_fft_scalar)(w*tonal->inmem[i]); |
| 257 | in[i].i = (kiss_fft_scalar)(w*tonal->inmem[N2+i]); |
| 258 | in[N-i-1].r = (kiss_fft_scalar)(w*tonal->inmem[N-i-1]); |
| 259 | in[N-i-1].i = (kiss_fft_scalar)(w*tonal->inmem[N+N2-i-1]); |
| 260 | } |
| 261 | OPUS_MOVE(tonal->inmem, tonal->inmem+ANALYSIS_BUF_SIZE-240, 240); |
| 262 | remaining = len - (ANALYSIS_BUF_SIZE-tonal->mem_fill); |
| 263 | downmix(x, &tonal->inmem[240], remaining, offset+ANALYSIS_BUF_SIZE-tonal->mem_fill, c1, c2, C); |
| 264 | tonal->mem_fill = 240 + remaining; |
| 265 | opus_fft(kfft, in, out); |
| 266 | |
| 267 | for (i=1;i<N2;i++) |
| 268 | { |
| 269 | float X1r, X2r, X1i, X2i; |
| 270 | float angle, d_angle, d2_angle; |
| 271 | float angle2, d_angle2, d2_angle2; |
| 272 | float mod1, mod2, avg_mod; |
| 273 | X1r = (float)out[i].r+out[N-i].r; |
| 274 | X1i = (float)out[i].i-out[N-i].i; |
| 275 | X2r = (float)out[i].i+out[N-i].i; |
| 276 | X2i = (float)out[N-i].r-out[i].r; |
| 277 | |
| 278 | angle = (float)(.5f/M_PI)*fast_atan2f(X1i, X1r); |
| 279 | d_angle = angle - A[i]; |
| 280 | d2_angle = d_angle - dA[i]; |
| 281 | |
| 282 | angle2 = (float)(.5f/M_PI)*fast_atan2f(X2i, X2r); |
| 283 | d_angle2 = angle2 - angle; |
| 284 | d2_angle2 = d_angle2 - d_angle; |
| 285 | |
| 286 | mod1 = d2_angle - (float)floor(.5+d2_angle); |
| 287 | noisiness[i] = ABS16(mod1); |
| 288 | mod1 *= mod1; |
| 289 | mod1 *= mod1; |
| 290 | |
| 291 | mod2 = d2_angle2 - (float)floor(.5+d2_angle2); |
| 292 | noisiness[i] += ABS16(mod2); |
| 293 | mod2 *= mod2; |
| 294 | mod2 *= mod2; |
| 295 | |
| 296 | avg_mod = .25f*(d2A[i]+2.f*mod1+mod2); |
| 297 | tonality[i] = 1.f/(1.f+40.f*16.f*pi4*avg_mod)-.015f; |
| 298 | |
| 299 | A[i] = angle2; |
| 300 | dA[i] = d_angle2; |
| 301 | d2A[i] = mod2; |
| 302 | } |
| 303 | |
| 304 | frame_tonality = 0; |
| 305 | max_frame_tonality = 0; |
| 306 | /*tw_sum = 0;*/ |
| 307 | info->activity = 0; |
| 308 | frame_noisiness = 0; |
| 309 | frame_stationarity = 0; |
| 310 | if (!tonal->count) |
| 311 | { |
| 312 | for (b=0;b<NB_TBANDS;b++) |
| 313 | { |
| 314 | tonal->lowE[b] = 1e10; |
| 315 | tonal->highE[b] = -1e10; |
| 316 | } |
| 317 | } |
| 318 | relativeE = 0; |
| 319 | frame_loudness = 0; |
| 320 | for (b=0;b<NB_TBANDS;b++) |
| 321 | { |
| 322 | float E=0, tE=0, nE=0; |
| 323 | float L1, L2; |
| 324 | float stationarity; |
| 325 | for (i=tbands[b];i<tbands[b+1];i++) |
| 326 | { |
| 327 | float binE = out[i].r*(float)out[i].r + out[N-i].r*(float)out[N-i].r |
| 328 | + out[i].i*(float)out[i].i + out[N-i].i*(float)out[N-i].i; |
| 329 | #ifdef FIXED_POINT |
| 330 | /* FIXME: It's probably best to change the BFCC filter initial state instead */ |
| 331 | binE *= 5.55e-17f; |
| 332 | #endif |
| 333 | E += binE; |
| 334 | tE += binE*tonality[i]; |
| 335 | nE += binE*2.f*(.5f-noisiness[i]); |
| 336 | } |
| 337 | tonal->E[tonal->E_count][b] = E; |
| 338 | frame_noisiness += nE/(1e-15f+E); |
| 339 | |
| 340 | frame_loudness += (float)sqrt(E+1e-10f); |
| 341 | logE[b] = (float)log(E+1e-10f); |
| 342 | tonal->lowE[b] = MIN32(logE[b], tonal->lowE[b]+.01f); |
| 343 | tonal->highE[b] = MAX32(logE[b], tonal->highE[b]-.1f); |
| 344 | if (tonal->highE[b] < tonal->lowE[b]+1.f) |
| 345 | { |
| 346 | tonal->highE[b]+=.5f; |
| 347 | tonal->lowE[b]-=.5f; |
| 348 | } |
| 349 | relativeE += (logE[b]-tonal->lowE[b])/(1e-15f+tonal->highE[b]-tonal->lowE[b]); |
| 350 | |
| 351 | L1=L2=0; |
| 352 | for (i=0;i<NB_FRAMES;i++) |
| 353 | { |
| 354 | L1 += (float)sqrt(tonal->E[i][b]); |
| 355 | L2 += tonal->E[i][b]; |
| 356 | } |
| 357 | |
| 358 | stationarity = MIN16(0.99f,L1/(float)sqrt(1e-15+NB_FRAMES*L2)); |
| 359 | stationarity *= stationarity; |
| 360 | stationarity *= stationarity; |
| 361 | frame_stationarity += stationarity; |
| 362 | /*band_tonality[b] = tE/(1e-15+E)*/; |
| 363 | band_tonality[b] = MAX16(tE/(1e-15f+E), stationarity*tonal->prev_band_tonality[b]); |
| 364 | #if 0 |
| 365 | if (b>=NB_TONAL_SKIP_BANDS) |
| 366 | { |
| 367 | frame_tonality += tweight[b]*band_tonality[b]; |
| 368 | tw_sum += tweight[b]; |
| 369 | } |
| 370 | #else |
| 371 | frame_tonality += band_tonality[b]; |
| 372 | if (b>=NB_TBANDS-NB_TONAL_SKIP_BANDS) |
| 373 | frame_tonality -= band_tonality[b-NB_TBANDS+NB_TONAL_SKIP_BANDS]; |
| 374 | #endif |
| 375 | max_frame_tonality = MAX16(max_frame_tonality, (1.f+.03f*(b-NB_TBANDS))*frame_tonality); |
| 376 | slope += band_tonality[b]*(b-8); |
| 377 | /*printf("%f %f ", band_tonality[b], stationarity);*/ |
| 378 | tonal->prev_band_tonality[b] = band_tonality[b]; |
| 379 | } |
| 380 | |
| 381 | bandwidth_mask = 0; |
| 382 | bandwidth = 0; |
| 383 | maxE = 0; |
| 384 | noise_floor = 5.7e-4f/(1<<(IMAX(0,lsb_depth-8))); |
| 385 | #ifdef FIXED_POINT |
| 386 | noise_floor *= 1<<(15+SIG_SHIFT); |
| 387 | #endif |
| 388 | noise_floor *= noise_floor; |
| 389 | for (b=0;b<NB_TOT_BANDS;b++) |
| 390 | { |
| 391 | float E=0; |
| 392 | int band_start, band_end; |
| 393 | /* Keep a margin of 300 Hz for aliasing */ |
| 394 | band_start = extra_bands[b]; |
| 395 | band_end = extra_bands[b+1]; |
| 396 | for (i=band_start;i<band_end;i++) |
| 397 | { |
| 398 | float binE = out[i].r*(float)out[i].r + out[N-i].r*(float)out[N-i].r |
| 399 | + out[i].i*(float)out[i].i + out[N-i].i*(float)out[N-i].i; |
| 400 | E += binE; |
| 401 | } |
| 402 | maxE = MAX32(maxE, E); |
| 403 | tonal->meanE[b] = MAX32((1-alphaE2)*tonal->meanE[b], E); |
| 404 | E = MAX32(E, tonal->meanE[b]); |
| 405 | /* Use a simple follower with 13 dB/Bark slope for spreading function */ |
| 406 | bandwidth_mask = MAX32(.05f*bandwidth_mask, E); |
| 407 | /* Consider the band "active" only if all these conditions are met: |
| 408 | 1) less than 10 dB below the simple follower |
| 409 | 2) less than 90 dB below the peak band (maximal masking possible considering |
| 410 | both the ATH and the loudness-dependent slope of the spreading function) |
| 411 | 3) above the PCM quantization noise floor |
| 412 | */ |
| 413 | if (E>.1*bandwidth_mask && E*1e9f > maxE && E > noise_floor*(band_end-band_start)) |
| 414 | bandwidth = b; |
| 415 | } |
| 416 | if (tonal->count<=2) |
| 417 | bandwidth = 20; |
| 418 | frame_loudness = 20*(float)log10(frame_loudness); |
| 419 | tonal->Etracker = MAX32(tonal->Etracker-.03f, frame_loudness); |
| 420 | tonal->lowECount *= (1-alphaE); |
| 421 | if (frame_loudness < tonal->Etracker-30) |
| 422 | tonal->lowECount += alphaE; |
| 423 | |
| 424 | for (i=0;i<8;i++) |
| 425 | { |
| 426 | float sum=0; |
| 427 | for (b=0;b<16;b++) |
| 428 | sum += dct_table[i*16+b]*logE[b]; |
| 429 | BFCC[i] = sum; |
| 430 | } |
| 431 | |
| 432 | frame_stationarity /= NB_TBANDS; |
| 433 | relativeE /= NB_TBANDS; |
| 434 | if (tonal->count<10) |
| 435 | relativeE = .5; |
| 436 | frame_noisiness /= NB_TBANDS; |
| 437 | #if 1 |
| 438 | info->activity = frame_noisiness + (1-frame_noisiness)*relativeE; |
| 439 | #else |
| 440 | info->activity = .5*(1+frame_noisiness-frame_stationarity); |
| 441 | #endif |
| 442 | frame_tonality = (max_frame_tonality/(NB_TBANDS-NB_TONAL_SKIP_BANDS)); |
| 443 | frame_tonality = MAX16(frame_tonality, tonal->prev_tonality*.8f); |
| 444 | tonal->prev_tonality = frame_tonality; |
| 445 | |
| 446 | slope /= 8*8; |
| 447 | info->tonality_slope = slope; |
| 448 | |
| 449 | tonal->E_count = (tonal->E_count+1)%NB_FRAMES; |
| 450 | tonal->count++; |
| 451 | info->tonality = frame_tonality; |
| 452 | |
| 453 | for (i=0;i<4;i++) |
| 454 | features[i] = -0.12299f*(BFCC[i]+tonal->mem[i+24]) + 0.49195f*(tonal->mem[i]+tonal->mem[i+16]) + 0.69693f*tonal->mem[i+8] - 1.4349f*tonal->cmean[i]; |
| 455 | |
| 456 | for (i=0;i<4;i++) |
| 457 | tonal->cmean[i] = (1-alpha)*tonal->cmean[i] + alpha*BFCC[i]; |
| 458 | |
| 459 | for (i=0;i<4;i++) |
| 460 | features[4+i] = 0.63246f*(BFCC[i]-tonal->mem[i+24]) + 0.31623f*(tonal->mem[i]-tonal->mem[i+16]); |
| 461 | for (i=0;i<3;i++) |
| 462 | features[8+i] = 0.53452f*(BFCC[i]+tonal->mem[i+24]) - 0.26726f*(tonal->mem[i]+tonal->mem[i+16]) -0.53452f*tonal->mem[i+8]; |
| 463 | |
| 464 | if (tonal->count > 5) |
| 465 | { |
| 466 | for (i=0;i<9;i++) |
| 467 | tonal->std[i] = (1-alpha)*tonal->std[i] + alpha*features[i]*features[i]; |
| 468 | } |
| 469 | |
| 470 | for (i=0;i<8;i++) |
| 471 | { |
| 472 | tonal->mem[i+24] = tonal->mem[i+16]; |
| 473 | tonal->mem[i+16] = tonal->mem[i+8]; |
| 474 | tonal->mem[i+8] = tonal->mem[i]; |
| 475 | tonal->mem[i] = BFCC[i]; |
| 476 | } |
| 477 | for (i=0;i<9;i++) |
| 478 | features[11+i] = (float)sqrt(tonal->std[i]); |
| 479 | features[20] = info->tonality; |
| 480 | features[21] = info->activity; |
| 481 | features[22] = frame_stationarity; |
| 482 | features[23] = info->tonality_slope; |
| 483 | features[24] = tonal->lowECount; |
| 484 | |
| 485 | #ifndef DISABLE_FLOAT_API |
| 486 | mlp_process(&net, features, frame_probs); |
| 487 | frame_probs[0] = .5f*(frame_probs[0]+1); |
| 488 | /* Curve fitting between the MLP probability and the actual probability */ |
| 489 | frame_probs[0] = .01f + 1.21f*frame_probs[0]*frame_probs[0] - .23f*(float)pow(frame_probs[0], 10); |
| 490 | /* Probability of active audio (as opposed to silence) */ |
| 491 | frame_probs[1] = .5f*frame_probs[1]+.5f; |
| 492 | /* Consider that silence has a 50-50 probability. */ |
| 493 | frame_probs[0] = frame_probs[1]*frame_probs[0] + (1-frame_probs[1])*.5f; |
| 494 | |
| 495 | /*printf("%f %f ", frame_probs[0], frame_probs[1]);*/ |
| 496 | { |
| 497 | /* Probability of state transition */ |
| 498 | float tau; |
| 499 | /* Represents independence of the MLP probabilities, where |
| 500 | beta=1 means fully independent. */ |
| 501 | float beta; |
| 502 | /* Denormalized probability of speech (p0) and music (p1) after update */ |
| 503 | float p0, p1; |
| 504 | /* Probabilities for "all speech" and "all music" */ |
| 505 | float s0, m0; |
| 506 | /* Probability sum for renormalisation */ |
| 507 | float psum; |
| 508 | /* Instantaneous probability of speech and music, with beta pre-applied. */ |
| 509 | float speech0; |
| 510 | float music0; |
| 511 | |
| 512 | /* One transition every 3 minutes of active audio */ |
| 513 | tau = .00005f*frame_probs[1]; |
| 514 | beta = .05f; |
| 515 | if (1) { |
| 516 | /* Adapt beta based on how "unexpected" the new prob is */ |
| 517 | float p, q; |
| 518 | p = MAX16(.05f,MIN16(.95f,frame_probs[0])); |
| 519 | q = MAX16(.05f,MIN16(.95f,tonal->music_prob)); |
| 520 | beta = .01f+.05f*ABS16(p-q)/(p*(1-q)+q*(1-p)); |
| 521 | } |
| 522 | /* p0 and p1 are the probabilities of speech and music at this frame |
| 523 | using only information from previous frame and applying the |
| 524 | state transition model */ |
| 525 | p0 = (1-tonal->music_prob)*(1-tau) + tonal->music_prob *tau; |
| 526 | p1 = tonal->music_prob *(1-tau) + (1-tonal->music_prob)*tau; |
| 527 | /* We apply the current probability with exponent beta to work around |
| 528 | the fact that the probability estimates aren't independent. */ |
| 529 | p0 *= (float)pow(1-frame_probs[0], beta); |
| 530 | p1 *= (float)pow(frame_probs[0], beta); |
| 531 | /* Normalise the probabilities to get the Marokv probability of music. */ |
| 532 | tonal->music_prob = p1/(p0+p1); |
| 533 | info->music_prob = tonal->music_prob; |
| 534 | |
| 535 | /* This chunk of code deals with delayed decision. */ |
| 536 | psum=1e-20f; |
| 537 | /* Instantaneous probability of speech and music, with beta pre-applied. */ |
| 538 | speech0 = (float)pow(1-frame_probs[0], beta); |
| 539 | music0 = (float)pow(frame_probs[0], beta); |
| 540 | if (tonal->count==1) |
| 541 | { |
| 542 | tonal->pspeech[0]=.5; |
| 543 | tonal->pmusic [0]=.5; |
| 544 | } |
| 545 | /* Updated probability of having only speech (s0) or only music (m0), |
| 546 | before considering the new observation. */ |
| 547 | s0 = tonal->pspeech[0] + tonal->pspeech[1]; |
| 548 | m0 = tonal->pmusic [0] + tonal->pmusic [1]; |
| 549 | /* Updates s0 and m0 with instantaneous probability. */ |
| 550 | tonal->pspeech[0] = s0*(1-tau)*speech0; |
| 551 | tonal->pmusic [0] = m0*(1-tau)*music0; |
| 552 | /* Propagate the transition probabilities */ |
| 553 | for (i=1;i<DETECT_SIZE-1;i++) |
| 554 | { |
| 555 | tonal->pspeech[i] = tonal->pspeech[i+1]*speech0; |
| 556 | tonal->pmusic [i] = tonal->pmusic [i+1]*music0; |
| 557 | } |
| 558 | /* Probability that the latest frame is speech, when all the previous ones were music. */ |
| 559 | tonal->pspeech[DETECT_SIZE-1] = m0*tau*speech0; |
| 560 | /* Probability that the latest frame is music, when all the previous ones were speech. */ |
| 561 | tonal->pmusic [DETECT_SIZE-1] = s0*tau*music0; |
| 562 | |
| 563 | /* Renormalise probabilities to 1 */ |
| 564 | for (i=0;i<DETECT_SIZE;i++) |
| 565 | psum += tonal->pspeech[i] + tonal->pmusic[i]; |
| 566 | psum = 1.f/psum; |
| 567 | for (i=0;i<DETECT_SIZE;i++) |
| 568 | { |
| 569 | tonal->pspeech[i] *= psum; |
| 570 | tonal->pmusic [i] *= psum; |
| 571 | } |
| 572 | psum = tonal->pmusic[0]; |
| 573 | for (i=1;i<DETECT_SIZE;i++) |
| 574 | psum += tonal->pspeech[i]; |
| 575 | |
| 576 | /* Estimate our confidence in the speech/music decisions */ |
| 577 | if (frame_probs[1]>.75) |
| 578 | { |
| 579 | if (tonal->music_prob>.9) |
| 580 | { |
| 581 | float adapt; |
| 582 | adapt = 1.f/(++tonal->music_confidence_count); |
| 583 | tonal->music_confidence_count = IMIN(tonal->music_confidence_count, 500); |
| 584 | tonal->music_confidence += adapt*MAX16(-.2f,frame_probs[0]-tonal->music_confidence); |
| 585 | } |
| 586 | if (tonal->music_prob<.1) |
| 587 | { |
| 588 | float adapt; |
| 589 | adapt = 1.f/(++tonal->speech_confidence_count); |
| 590 | tonal->speech_confidence_count = IMIN(tonal->speech_confidence_count, 500); |
| 591 | tonal->speech_confidence += adapt*MIN16(.2f,frame_probs[0]-tonal->speech_confidence); |
| 592 | } |
| 593 | } else { |
| 594 | if (tonal->music_confidence_count==0) |
| 595 | tonal->music_confidence = .9f; |
| 596 | if (tonal->speech_confidence_count==0) |
| 597 | tonal->speech_confidence = .1f; |
| 598 | } |
| 599 | } |
| 600 | if (tonal->last_music != (tonal->music_prob>.5f)) |
| 601 | tonal->last_transition=0; |
| 602 | tonal->last_music = tonal->music_prob>.5f; |
| 603 | #else |
| 604 | info->music_prob = 0; |
| 605 | #endif |
| 606 | /*for (i=0;i<25;i++) |
| 607 | printf("%f ", features[i]); |
| 608 | printf("\n");*/ |
| 609 | |
| 610 | info->bandwidth = bandwidth; |
| 611 | /*printf("%d %d\n", info->bandwidth, info->opus_bandwidth);*/ |
| 612 | info->noisiness = frame_noisiness; |
| 613 | info->valid = 1; |
| 614 | if (info_out!=NULL) |
| 615 | OPUS_COPY(info_out, info, 1); |
| 616 | RESTORE_STACK; |
| 617 | } |
| 618 | |
| 619 | void run_analysis(TonalityAnalysisState *analysis, const CELTMode *celt_mode, const void *analysis_pcm, |
| 620 | int analysis_frame_size, int frame_size, int c1, int c2, int C, opus_int32 Fs, |
| 621 | int lsb_depth, downmix_func downmix, AnalysisInfo *analysis_info) |
| 622 | { |
| 623 | int offset; |
| 624 | int pcm_len; |
| 625 | |
| 626 | if (analysis_pcm != NULL) |
| 627 | { |
| 628 | /* Avoid overflow/wrap-around of the analysis buffer */ |
| 629 | analysis_frame_size = IMIN((DETECT_SIZE-5)*Fs/100, analysis_frame_size); |
| 630 | |
| 631 | pcm_len = analysis_frame_size - analysis->analysis_offset; |
| 632 | offset = analysis->analysis_offset; |
| 633 | do { |
| 634 | tonality_analysis(analysis, NULL, celt_mode, analysis_pcm, IMIN(480, pcm_len), offset, c1, c2, C, lsb_depth, downmix); |
| 635 | offset += 480; |
| 636 | pcm_len -= 480; |
| 637 | } while (pcm_len>0); |
| 638 | analysis->analysis_offset = analysis_frame_size; |
| 639 | |
| 640 | analysis->analysis_offset -= frame_size; |
| 641 | } |
| 642 | |
| 643 | analysis_info->valid = 0; |
| 644 | tonality_get_info(analysis, analysis_info, frame_size); |
| 645 | } |