diff options
Diffstat (limited to 'src/modules/audio_processing/ns/main/source/nsx_core.c')
| -rw-r--r-- | src/modules/audio_processing/ns/main/source/nsx_core.c | 2493 |
1 files changed, 0 insertions, 2493 deletions
diff --git a/src/modules/audio_processing/ns/main/source/nsx_core.c b/src/modules/audio_processing/ns/main/source/nsx_core.c deleted file mode 100644 index 01d3e54080..0000000000 --- a/src/modules/audio_processing/ns/main/source/nsx_core.c +++ /dev/null @@ -1,2493 +0,0 @@ -/* - * Copyright (c) 2011 The WebRTC project authors. All Rights Reserved. - * - * Use of this source code is governed by a BSD-style license - * that can be found in the LICENSE file in the root of the source - * tree. An additional intellectual property rights grant can be found - * in the file PATENTS. All contributing project authors may - * be found in the AUTHORS file in the root of the source tree. - */ - -#include "noise_suppression_x.h" - -#include <assert.h> -#include <math.h> -#include <string.h> -#include <stdlib.h> - -#include "nsx_core.h" - -// Skip first frequency bins during estimation. (0 <= value < 64) -static const int kStartBand = 5; - -// Rounding -static const WebRtc_Word16 kRoundTable[16] = {0, 1, 2, 4, 8, 16, 32, 64, 128, 256, 512, 1024, - 2048, 4096, 8192, 16384}; - -// Constants to compensate for shifting signal log(2^shifts). -static const WebRtc_Word16 kLogTable[9] = {0, 177, 355, 532, 710, 887, 1065, 1242, 1420}; - -static const WebRtc_Word16 kCounterDiv[201] = {32767, 16384, 10923, 8192, 6554, 5461, 4681, - 4096, 3641, 3277, 2979, 2731, 2521, 2341, 2185, 2048, 1928, 1820, 1725, 1638, 1560, - 1489, 1425, 1365, 1311, 1260, 1214, 1170, 1130, 1092, 1057, 1024, 993, 964, 936, 910, - 886, 862, 840, 819, 799, 780, 762, 745, 728, 712, 697, 683, 669, 655, 643, 630, 618, - 607, 596, 585, 575, 565, 555, 546, 537, 529, 520, 512, 504, 496, 489, 482, 475, 468, - 462, 455, 449, 443, 437, 431, 426, 420, 415, 410, 405, 400, 395, 390, 386, 381, 377, - 372, 368, 364, 360, 356, 352, 349, 345, 341, 338, 334, 331, 328, 324, 321, 318, 315, - 312, 309, 306, 303, 301, 298, 295, 293, 290, 287, 285, 282, 280, 278, 275, 273, 271, - 269, 266, 264, 262, 260, 258, 256, 254, 252, 250, 248, 246, 245, 243, 241, 239, 237, - 236, 234, 232, 231, 229, 228, 226, 224, 223, 221, 220, 218, 217, 216, 214, 213, 211, - 210, 209, 207, 206, 205, 204, 202, 201, 200, 199, 197, 196, 195, 194, 193, 192, 191, - 189, 188, 187, 186, 185, 184, 183, 182, 181, 180, 179, 178, 177, 176, 175, 174, 173, - 172, 172, 171, 170, 169, 168, 167, 166, 165, 165, 164, 163}; - -static const WebRtc_Word16 kLogTableFrac[256] = { - 0, 1, 3, 4, 6, 7, 9, 10, 11, 13, 14, 16, 17, 18, 20, 21, - 22, 24, 25, 26, 28, 29, 30, 32, 33, 34, 36, 37, 38, 40, 41, 42, - 44, 45, 46, 47, 49, 50, 51, 52, 54, 55, 56, 57, 59, 60, 61, 62, - 63, 65, 66, 67, 68, 69, 71, 72, 73, 74, 75, 77, 78, 79, 80, 81, - 82, 84, 85, 86, 87, 88, 89, 90, 92, 93, 94, 95, 96, 97, 98, 99, - 100, 102, 103, 104, 105, 106, 107, 108, 109, 110, 111, 112, 113, 114, 116, 117, - 118, 119, 120, 121, 122, 123, 124, 125, 126, 127, 128, 129, 130, 131, 132, 133, - 134, 135, 136, 137, 138, 139, 140, 141, 142, 143, 144, 145, 146, 147, 148, 149, - 150, 151, 152, 153, 154, 155, 155, 156, 157, 158, 159, 160, 161, 162, 163, 164, - 165, 166, 167, 168, 169, 169, 170, 171, 172, 173, 174, 175, 176, 177, 178, 178, - 179, 180, 181, 182, 183, 184, 185, 185, 186, 187, 188, 189, 190, 191, 192, 192, - 193, 194, 195, 196, 197, 198, 198, 199, 200, 201, 202, 203, 203, 204, 205, 206, - 207, 208, 208, 209, 210, 211, 212, 212, 213, 214, 215, 216, 216, 217, 218, 219, - 220, 220, 221, 222, 223, 224, 224, 225, 226, 227, 228, 228, 229, 230, 231, 231, - 232, 233, 234, 234, 235, 236, 237, 238, 238, 239, 240, 241, 241, 242, 243, 244, - 244, 245, 246, 247, 247, 248, 249, 249, 250, 251, 252, 252, 253, 254, 255, 255 -}; - -static const WebRtc_Word16 kPowTableFrac[1024] = { - 0, 1, 1, 2, 3, 3, 4, 5, - 6, 6, 7, 8, 8, 9, 10, 10, - 11, 12, 13, 13, 14, 15, 15, 16, - 17, 17, 18, 19, 20, 20, 21, 22, - 22, 23, 24, 25, 25, 26, 27, 27, - 28, 29, 30, 30, 31, 32, 32, 33, - 34, 35, 35, 36, 37, 37, 38, 39, - 40, 40, 41, 42, 42, 43, 44, 45, - 45, 46, 47, 48, 48, 49, 50, 50, - 51, 52, 53, 53, 54, 55, 56, 56, - 57, 58, 58, 59, 60, 61, 61, 62, - 63, 64, 64, 65, 66, 67, 67, 68, - 69, 69, 70, 71, 72, 72, 73, 74, - 75, 75, 76, 77, 78, 78, 79, 80, - 81, 81, 82, 83, 84, 84, 85, 86, - 87, 87, 88, 89, 90, 90, 91, 92, - 93, 93, 94, 95, 96, 96, 97, 98, - 99, 100, 100, 101, 102, 103, 103, 104, - 105, 106, 106, 107, 108, 109, 109, 110, - 111, 112, 113, 113, 114, 115, 116, 116, - 117, 118, 119, 119, 120, 121, 122, 123, - 123, 124, 125, 126, 126, 127, 128, 129, - 130, 130, 131, 132, 133, 133, 134, 135, - 136, 137, 137, 138, 139, 140, 141, 141, - 142, 143, 144, 144, 145, 146, 147, 148, - 148, 149, 150, 151, 152, 152, 153, 154, - 155, 156, 156, 157, 158, 159, 160, 160, - 161, 162, 163, 164, 164, 165, 166, 167, - 168, 168, 169, 170, 171, 172, 173, 173, - 174, 175, 176, 177, 177, 178, 179, 180, - 181, 181, 182, 183, 184, 185, 186, 186, - 187, 188, 189, 190, 190, 191, 192, 193, - 194, 195, 195, 196, 197, 198, 199, 200, - 200, 201, 202, 203, 204, 205, 205, 206, - 207, 208, 209, 210, 210, 211, 212, 213, - 214, 215, 215, 216, 217, 218, 219, 220, - 220, 221, 222, 223, 224, 225, 225, 226, - 227, 228, 229, 230, 231, 231, 232, 233, - 234, 235, 236, 237, 237, 238, 239, 240, - 241, 242, 243, 243, 244, 245, 246, 247, - 248, 249, 249, 250, 251, 252, 253, 254, - 255, 255, 256, 257, 258, 259, 260, 261, - 262, 262, 263, 264, 265, 266, 267, 268, - 268, 269, 270, 271, 272, 273, 274, 275, - 276, 276, 277, 278, 279, 280, 281, 282, - 283, 283, 284, 285, 286, 287, 288, 289, - 290, 291, 291, 292, 293, 294, 295, 296, - 297, 298, 299, 299, 300, 301, 302, 303, - 304, 305, 306, 307, 308, 308, 309, 310, - 311, 312, 313, 314, 315, 316, 317, 318, - 318, 319, 320, 321, 322, 323, 324, 325, - 326, 327, 328, 328, 329, 330, 331, 332, - 333, 334, 335, 336, 337, 338, 339, 339, - 340, 341, 342, 343, 344, 345, 346, 347, - 348, 349, 350, 351, 352, 352, 353, 354, - 355, 356, 357, 358, 359, 360, 361, 362, - 363, 364, 365, 366, 367, 367, 368, 369, - 370, 371, 372, 373, 374, 375, 376, 377, - 378, 379, 380, 381, 382, 383, 384, 385, - 385, 386, 387, 388, 389, 390, 391, 392, - 393, 394, 395, 396, 397, 398, 399, 400, - 401, 402, 403, 404, 405, 406, 407, 408, - 409, 410, 410, 411, 412, 413, 414, 415, - 416, 417, 418, 419, 420, 421, 422, 423, - 424, 425, 426, 427, 428, 429, 430, 431, - 432, 433, 434, 435, 436, 437, 438, 439, - 440, 441, 442, 443, 444, 445, 446, 447, - 448, 449, 450, 451, 452, 453, 454, 455, - 456, 457, 458, 459, 460, 461, 462, 463, - 464, 465, 466, 467, 468, 469, 470, 471, - 472, 473, 474, 475, 476, 477, 478, 479, - 480, 481, 482, 483, 484, 485, 486, 487, - 488, 489, 490, 491, 492, 493, 494, 495, - 496, 498, 499, 500, 501, 502, 503, 504, - 505, 506, 507, 508, 509, 510, 511, 512, - 513, 514, 515, 516, 517, 518, 519, 520, - 521, 522, 523, 525, 526, 527, 528, 529, - 530, 531, 532, 533, 534, 535, 536, 537, - 538, 539, 540, 541, 542, 544, 545, 546, - 547, 548, 549, 550, 551, 552, 553, 554, - 555, 556, 557, 558, 560, 561, 562, 563, - 564, 565, 566, 567, 568, 569, 570, 571, - 572, 574, 575, 576, 577, 578, 579, 580, - 581, 582, 583, 584, 585, 587, 588, 589, - 590, 591, 592, 593, 594, 595, 596, 597, - 599, 600, 601, 602, 603, 604, 605, 606, - 607, 608, 610, 611, 612, 613, 614, 615, - 616, 617, 618, 620, 621, 622, 623, 624, - 625, 626, 627, 628, 630, 631, 632, 633, - 634, 635, 636, 637, 639, 640, 641, 642, - 643, 644, 645, 646, 648, 649, 650, 651, - 652, 653, 654, 656, 657, 658, 659, 660, - 661, 662, 664, 665, 666, 667, 668, 669, - 670, 672, 673, 674, 675, 676, 677, 678, - 680, 681, 682, 683, 684, 685, 687, 688, - 689, 690, 691, 692, 693, 695, 696, 697, - 698, 699, 700, 702, 703, 704, 705, 706, - 708, 709, 710, 711, 712, 713, 715, 716, - 717, 718, 719, 720, 722, 723, 724, 725, - 726, 728, 729, 730, 731, 732, 733, 735, - 736, 737, 738, 739, 741, 742, 743, 744, - 745, 747, 748, 749, 750, 751, 753, 754, - 755, 756, 757, 759, 760, 761, 762, 763, - 765, 766, 767, 768, 770, 771, 772, 773, - 774, 776, 777, 778, 779, 780, 782, 783, - 784, 785, 787, 788, 789, 790, 792, 793, - 794, 795, 796, 798, 799, 800, 801, 803, - 804, 805, 806, 808, 809, 810, 811, 813, - 814, 815, 816, 818, 819, 820, 821, 823, - 824, 825, 826, 828, 829, 830, 831, 833, - 834, 835, 836, 838, 839, 840, 841, 843, - 844, 845, 846, 848, 849, 850, 851, 853, - 854, 855, 857, 858, 859, 860, 862, 863, - 864, 866, 867, 868, 869, 871, 872, 873, - 874, 876, 877, 878, 880, 881, 882, 883, - 885, 886, 887, 889, 890, 891, 893, 894, - 895, 896, 898, 899, 900, 902, 903, 904, - 906, 907, 908, 909, 911, 912, 913, 915, - 916, 917, 919, 920, 921, 923, 924, 925, - 927, 928, 929, 931, 932, 933, 935, 936, - 937, 938, 940, 941, 942, 944, 945, 946, - 948, 949, 950, 952, 953, 955, 956, 957, - 959, 960, 961, 963, 964, 965, 967, 968, - 969, 971, 972, 973, 975, 976, 977, 979, - 980, 981, 983, 984, 986, 987, 988, 990, - 991, 992, 994, 995, 996, 998, 999, 1001, - 1002, 1003, 1005, 1006, 1007, 1009, 1010, 1012, - 1013, 1014, 1016, 1017, 1018, 1020, 1021, 1023 -}; - -static const WebRtc_Word16 kIndicatorTable[17] = {0, 2017, 3809, 5227, 6258, 6963, 7424, 7718, - 7901, 8014, 8084, 8126, 8152, 8168, 8177, 8183, 8187}; - -// hybrib Hanning & flat window -static const WebRtc_Word16 kBlocks80w128x[128] = { - 0, 536, 1072, 1606, 2139, 2669, 3196, 3720, 4240, 4756, 5266, - 5771, 6270, 6762, 7246, 7723, 8192, 8652, 9102, 9543, 9974, 10394, - 10803, 11200, 11585, 11958, 12318, 12665, 12998, 13318, 13623, 13913, 14189, - 14449, 14694, 14924, 15137, 15334, 15515, 15679, 15826, 15956, 16069, 16165, - 16244, 16305, 16349, 16375, 16384, 16384, 16384, 16384, 16384, 16384, 16384, - 16384, 16384, 16384, 16384, 16384, 16384, 16384, 16384, 16384, 16384, 16384, - 16384, 16384, 16384, 16384, 16384, 16384, 16384, 16384, 16384, 16384, 16384, - 16384, 16384, 16384, 16384, 16375, 16349, 16305, 16244, 16165, 16069, 15956, - 15826, 15679, 15515, 15334, 15137, 14924, 14694, 14449, 14189, 13913, 13623, - 13318, 12998, 12665, 12318, 11958, 11585, 11200, 10803, 10394, 9974, 9543, - 9102, 8652, 8192, 7723, 7246, 6762, 6270, 5771, 5266, 4756, 4240, - 3720, 3196, 2669, 2139, 1606, 1072, 536 -}; - -// hybrib Hanning & flat window -static const WebRtc_Word16 kBlocks160w256x[256] = { - 0, 268, 536, 804, 1072, 1339, 1606, 1872, - 2139, 2404, 2669, 2933, 3196, 3459, 3720, 3981, - 4240, 4499, 4756, 5012, 5266, 5520, 5771, 6021, - 6270, 6517, 6762, 7005, 7246, 7486, 7723, 7959, - 8192, 8423, 8652, 8878, 9102, 9324, 9543, 9760, - 9974, 10185, 10394, 10600, 10803, 11003, 11200, 11394, -11585, 11773, 11958, 12140, 12318, 12493, 12665, 12833, -12998, 13160, 13318, 13472, 13623, 13770, 13913, 14053, -14189, 14321, 14449, 14574, 14694, 14811, 14924, 15032, -15137, 15237, 15334, 15426, 15515, 15599, 15679, 15754, -15826, 15893, 15956, 16015, 16069, 16119, 16165, 16207, -16244, 16277, 16305, 16329, 16349, 16364, 16375, 16382, -16384, 16384, 16384, 16384, 16384, 16384, 16384, 16384, -16384, 16384, 16384, 16384, 16384, 16384, 16384, 16384, -16384, 16384, 16384, 16384, 16384, 16384, 16384, 16384, -16384, 16384, 16384, 16384, 16384, 16384, 16384, 16384, -16384, 16384, 16384, 16384, 16384, 16384, 16384, 16384, -16384, 16384, 16384, 16384, 16384, 16384, 16384, 16384, -16384, 16384, 16384, 16384, 16384, 16384, 16384, 16384, -16384, 16384, 16384, 16384, 16384, 16384, 16384, 16384, -16384, 16382, 16375, 16364, 16349, 16329, 16305, 16277, -16244, 16207, 16165, 16119, 16069, 16015, 15956, 15893, -15826, 15754, 15679, 15599, 15515, 15426, 15334, 15237, -15137, 15032, 14924, 14811, 14694, 14574, 14449, 14321, -14189, 14053, 13913, 13770, 13623, 13472, 13318, 13160, -12998, 12833, 12665, 12493, 12318, 12140, 11958, 11773, -11585, 11394, 11200, 11003, 10803, 10600, 10394, 10185, - 9974, 9760, 9543, 9324, 9102, 8878, 8652, 8423, - 8192, 7959, 7723, 7486, 7246, 7005, 6762, 6517, - 6270, 6021, 5771, 5520, 5266, 5012, 4756, 4499, - 4240, 3981, 3720, 3459, 3196, 2933, 2669, 2404, - 2139, 1872, 1606, 1339, 1072, 804, 536, 268 -}; - -// Gain factor table: Input value in Q8 and output value in Q13 -static const WebRtc_Word16 kFactor1Table[257] = { - 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8233, 8274, 8315, 8355, 8396, 8436, 8475, 8515, 8554, 8592, 8631, 8669, - 8707, 8745, 8783, 8820, 8857, 8894, 8931, 8967, 9003, 9039, 9075, 9111, 9146, 9181, - 9216, 9251, 9286, 9320, 9354, 9388, 9422, 9456, 9489, 9523, 9556, 9589, 9622, 9655, - 9687, 9719, 9752, 9784, 9816, 9848, 9879, 9911, 9942, 9973, 10004, 10035, 10066, - 10097, 10128, 10158, 10188, 10218, 10249, 10279, 10308, 10338, 10368, 10397, 10426, - 10456, 10485, 10514, 10543, 10572, 10600, 10629, 10657, 10686, 10714, 10742, 10770, - 10798, 10826, 10854, 10882, 10847, 10810, 10774, 10737, 10701, 10666, 10631, 10596, - 10562, 10527, 10494, 10460, 10427, 10394, 10362, 10329, 10297, 10266, 10235, 10203, - 10173, 10142, 10112, 10082, 10052, 10023, 9994, 9965, 9936, 9908, 9879, 9851, 9824, - 9796, 9769, 9742, 9715, 9689, 9662, 9636, 9610, 9584, 9559, 9534, 9508, 9484, 9459, - 9434, 9410, 9386, 9362, 9338, 9314, 9291, 9268, 9245, 9222, 9199, 9176, 9154, 9132, - 9110, 9088, 9066, 9044, 9023, 9002, 8980, 8959, 8939, 8918, 8897, 8877, 8857, 8836, - 8816, 8796, 8777, 8757, 8738, 8718, 8699, 8680, 8661, 8642, 8623, 8605, 8586, 8568, - 8550, 8532, 8514, 8496, 8478, 8460, 8443, 8425, 8408, 8391, 8373, 8356, 8339, 8323, - 8306, 8289, 8273, 8256, 8240, 8224, 8208, 8192 -}; - -// Gain factor table: Input value in Q8 and output value in Q13 -static const WebRtc_Word16 kFactor2Aggressiveness1[257] = { - 7577, 7577, 7577, 7577, 7577, 7577, - 7577, 7577, 7577, 7577, 7577, 7577, 7577, 7577, 7577, 7577, 7577, 7596, 7614, 7632, - 7650, 7667, 7683, 7699, 7715, 7731, 7746, 7761, 7775, 7790, 7804, 7818, 7832, 7845, - 7858, 7871, 7884, 7897, 7910, 7922, 7934, 7946, 7958, 7970, 7982, 7993, 8004, 8016, - 8027, 8038, 8049, 8060, 8070, 8081, 8091, 8102, 8112, 8122, 8132, 8143, 8152, 8162, - 8172, 8182, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192 -}; - -// Gain factor table: Input value in Q8 and output value in Q13 -static const WebRtc_Word16 kFactor2Aggressiveness2[257] = { - 7270, 7270, 7270, 7270, 7270, 7306, - 7339, 7369, 7397, 7424, 7448, 7472, 7495, 7517, 7537, 7558, 7577, 7596, 7614, 7632, - 7650, 7667, 7683, 7699, 7715, 7731, 7746, 7761, 7775, 7790, 7804, 7818, 7832, 7845, - 7858, 7871, 7884, 7897, 7910, 7922, 7934, 7946, 7958, 7970, 7982, 7993, 8004, 8016, - 8027, 8038, 8049, 8060, 8070, 8081, 8091, 8102, 8112, 8122, 8132, 8143, 8152, 8162, - 8172, 8182, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192 -}; - -// Gain factor table: Input value in Q8 and output value in Q13 -static const WebRtc_Word16 kFactor2Aggressiveness3[257] = { - 7184, 7184, 7184, 7229, 7270, 7306, - 7339, 7369, 7397, 7424, 7448, 7472, 7495, 7517, 7537, 7558, 7577, 7596, 7614, 7632, - 7650, 7667, 7683, 7699, 7715, 7731, 7746, 7761, 7775, 7790, 7804, 7818, 7832, 7845, - 7858, 7871, 7884, 7897, 7910, 7922, 7934, 7946, 7958, 7970, 7982, 7993, 8004, 8016, - 8027, 8038, 8049, 8060, 8070, 8081, 8091, 8102, 8112, 8122, 8132, 8143, 8152, 8162, - 8172, 8182, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, - 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192, 8192 -}; - -// sum of log2(i) from table index to inst->anaLen2 in Q5 -// Note that the first table value is invalid, since log2(0) = -infinity -static const WebRtc_Word16 kSumLogIndex[66] = { - 0, 22917, 22917, 22885, 22834, 22770, 22696, 22613, - 22524, 22428, 22326, 22220, 22109, 21994, 21876, 21754, - 21629, 21501, 21370, 21237, 21101, 20963, 20822, 20679, - 20535, 20388, 20239, 20089, 19937, 19783, 19628, 19470, - 19312, 19152, 18991, 18828, 18664, 18498, 18331, 18164, - 17994, 17824, 17653, 17480, 17306, 17132, 16956, 16779, - 16602, 16423, 16243, 16063, 15881, 15699, 15515, 15331, - 15146, 14960, 14774, 14586, 14398, 14209, 14019, 13829, - 13637, 13445 -}; - -// sum of log2(i)^2 from table index to inst->anaLen2 in Q2 -// Note that the first table value is invalid, since log2(0) = -infinity -static const WebRtc_Word16 kSumSquareLogIndex[66] = { - 0, 16959, 16959, 16955, 16945, 16929, 16908, 16881, - 16850, 16814, 16773, 16729, 16681, 16630, 16575, 16517, - 16456, 16392, 16325, 16256, 16184, 16109, 16032, 15952, - 15870, 15786, 15700, 15612, 15521, 15429, 15334, 15238, - 15140, 15040, 14938, 14834, 14729, 14622, 14514, 14404, - 14292, 14179, 14064, 13947, 13830, 13710, 13590, 13468, - 13344, 13220, 13094, 12966, 12837, 12707, 12576, 12444, - 12310, 12175, 12039, 11902, 11763, 11624, 11483, 11341, - 11198, 11054 -}; - -// log2(table index) in Q12 -// Note that the first table value is invalid, since log2(0) = -infinity -static const WebRtc_Word16 kLogIndex[129] = { - 0, 0, 4096, 6492, 8192, 9511, 10588, 11499, - 12288, 12984, 13607, 14170, 14684, 15157, 15595, 16003, - 16384, 16742, 17080, 17400, 17703, 17991, 18266, 18529, - 18780, 19021, 19253, 19476, 19691, 19898, 20099, 20292, - 20480, 20662, 20838, 21010, 21176, 21338, 21496, 21649, - 21799, 21945, 22087, 22226, 22362, 22495, 22625, 22752, - 22876, 22998, 23117, 23234, 23349, 23462, 23572, 23680, - 23787, 23892, 23994, 24095, 24195, 24292, 24388, 24483, - 24576, 24668, 24758, 24847, 24934, 25021, 25106, 25189, - 25272, 25354, 25434, 25513, 25592, 25669, 25745, 25820, - 25895, 25968, 26041, 26112, 26183, 26253, 26322, 26390, - 26458, 26525, 26591, 26656, 26721, 26784, 26848, 26910, - 26972, 27033, 27094, 27154, 27213, 27272, 27330, 27388, - 27445, 27502, 27558, 27613, 27668, 27722, 27776, 27830, - 27883, 27935, 27988, 28039, 28090, 28141, 28191, 28241, - 28291, 28340, 28388, 28437, 28484, 28532, 28579, 28626, - 28672 -}; - -// determinant of estimation matrix in Q0 corresponding to the log2 tables above -// Note that the first table value is invalid, since log2(0) = -infinity -static const WebRtc_Word16 kDeterminantEstMatrix[66] = { - 0, 29814, 25574, 22640, 20351, 18469, 16873, 15491, - 14277, 13199, 12233, 11362, 10571, 9851, 9192, 8587, - 8030, 7515, 7038, 6596, 6186, 5804, 5448, 5115, - 4805, 4514, 4242, 3988, 3749, 3524, 3314, 3116, - 2930, 2755, 2590, 2435, 2289, 2152, 2022, 1900, - 1785, 1677, 1575, 1478, 1388, 1302, 1221, 1145, - 1073, 1005, 942, 881, 825, 771, 721, 674, - 629, 587, 547, 510, 475, 442, 411, 382, - 355, 330 -}; - -void WebRtcNsx_UpdateNoiseEstimate(NsxInst_t *inst, int offset) -{ - WebRtc_Word32 tmp32no1 = 0; - WebRtc_Word32 tmp32no2 = 0; - - WebRtc_Word16 tmp16no1 = 0; - WebRtc_Word16 tmp16no2 = 0; - WebRtc_Word16 exp2Const = 11819; // Q13 - - int i = 0; - - tmp16no2 = WebRtcSpl_MaxValueW16(inst->noiseEstLogQuantile + offset, inst->magnLen); - inst->qNoise = 14 - - (int)WEBRTC_SPL_MUL_16_16_RSFT_WITH_ROUND(exp2Const, tmp16no2, 21); - for (i = 0; i < inst->magnLen; i++) - { - // inst->quantile[i]=exp(inst->lquantile[offset+i]); - // in Q21 - tmp32no2 = WEBRTC_SPL_MUL_16_16(exp2Const, inst->noiseEstLogQuantile[offset + i]); - tmp32no1 = (0x00200000 | (tmp32no2 & 0x001FFFFF)); - tmp16no1 = -(WebRtc_Word16)WEBRTC_SPL_RSHIFT_W32(tmp32no2, 21); - tmp16no1 += 21;// shift 21 to get result in Q0 - tmp16no1 -= (WebRtc_Word16)inst->qNoise; //shift to get result in Q(qNoise) - if (tmp16no1 > 0) - { - inst->noiseEstQuantile[i] = (WebRtc_Word16)WEBRTC_SPL_RSHIFT_W32(tmp32no1 + - kRoundTable[tmp16no1], tmp16no1); - } - else - { - inst->noiseEstQuantile[i] = (WebRtc_Word16)WEBRTC_SPL_LSHIFT_W32(tmp32no1, - -tmp16no1); - } - } -} - -void WebRtcNsx_CalcParametricNoiseEstimate(NsxInst_t *inst, - WebRtc_Word16 pink_noise_exp_avg, - WebRtc_Word32 pink_noise_num_avg, - int freq_index, - WebRtc_UWord32 *noise_estimate, - WebRtc_UWord32 *noise_estimate_avg) -{ - WebRtc_Word32 tmp32no1 = 0; - WebRtc_Word32 tmp32no2 = 0; - - WebRtc_Word16 int_part = 0; - WebRtc_Word16 frac_part = 0; - - // Use pink noise estimate - // noise_estimate = 2^(pinkNoiseNumerator + pinkNoiseExp * log2(j)) - assert(freq_index > 0); - tmp32no2 = WEBRTC_SPL_MUL_16_16(pink_noise_exp_avg, kLogIndex[freq_index]); // Q26 - tmp32no2 = WEBRTC_SPL_RSHIFT_W32(tmp32no2, 15); // Q11 - tmp32no1 = pink_noise_num_avg - tmp32no2; // Q11 - - // Calculate output: 2^tmp32no1 - // Output in Q(minNorm-stages) - tmp32no1 += WEBRTC_SPL_LSHIFT_W32((WebRtc_Word32)(inst->minNorm - inst->stages), 11); - if (tmp32no1 > 0) - { - int_part = (WebRtc_Word16)WEBRTC_SPL_RSHIFT_W32(tmp32no1, 11); - frac_part = (WebRtc_Word16)(tmp32no1 & 0x000007ff); // Q11 - // Piecewise linear approximation of 'b' in - // 2^(int_part+frac_part) = 2^int_part * (1 + b) - // 'b' is given in Q11 and below stored in frac_part. - if (WEBRTC_SPL_RSHIFT_W32(frac_part, 10)) - { - // Upper fractional part - tmp32no2 = WEBRTC_SPL_MUL_32_16(2048 - frac_part, 1244); // Q21 - tmp32no2 = 2048 - WEBRTC_SPL_RSHIFT_W32(tmp32no2, 10); - } - else - { - // Lower fractional part - tmp32no2 = WEBRTC_SPL_RSHIFT_W32(WEBRTC_SPL_MUL_32_16(frac_part, 804), 10); - } - // Shift fractional part to Q(minNorm-stages) - tmp32no2 = WEBRTC_SPL_SHIFT_W32(tmp32no2, int_part - 11); - *noise_estimate_avg = WEBRTC_SPL_LSHIFT_U32(1, int_part) + (WebRtc_UWord32)tmp32no2; - // Scale up to initMagnEst, which is not block averaged - *noise_estimate = (*noise_estimate_avg) * (WebRtc_UWord32)(inst->blockIndex + 1); - } -} - -// Initialize state -WebRtc_Word32 WebRtcNsx_InitCore(NsxInst_t *inst, WebRtc_UWord32 fs) -{ - int i; - - //check for valid pointer - if (inst == NULL) - { - return -1; - } - // - - // Initialization of struct - if (fs == 8000 || fs == 16000 || fs == 32000) - { - inst->fs = fs; - } else - { - return -1; - } - - if (fs == 8000) - { - inst->blockLen10ms = 80; - inst->anaLen = 128; - inst->stages = 7; - inst->window = kBlocks80w128x; - inst->thresholdLogLrt = 131072; //default threshold for LRT feature - inst->maxLrt = 0x0040000; - inst->minLrt = 52429; - } else if (fs == 16000) - { - inst->blockLen10ms = 160; - inst->anaLen = 256; - inst->stages = 8; - inst->window = kBlocks160w256x; - inst->thresholdLogLrt = 212644; //default threshold for LRT feature - inst->maxLrt = 0x0080000; - inst->minLrt = 104858; - } else if (fs == 32000) - { - inst->blockLen10ms = 160; - inst->anaLen = 256; - inst->stages = 8; - inst->window = kBlocks160w256x; - inst->thresholdLogLrt = 212644; //default threshold for LRT feature - inst->maxLrt = 0x0080000; - inst->minLrt = 104858; - } - inst->anaLen2 = WEBRTC_SPL_RSHIFT_W16(inst->anaLen, 1); - inst->magnLen = inst->anaLen2 + 1; - - WebRtcSpl_ZerosArrayW16(inst->analysisBuffer, ANAL_BLOCKL_MAX); - WebRtcSpl_ZerosArrayW16(inst->synthesisBuffer, ANAL_BLOCKL_MAX); - - // for HB processing - WebRtcSpl_ZerosArrayW16(inst->dataBufHBFX, ANAL_BLOCKL_MAX); - // for quantile noise estimation - WebRtcSpl_ZerosArrayW16(inst->noiseEstQuantile, HALF_ANAL_BLOCKL); - for (i = 0; i < SIMULT * HALF_ANAL_BLOCKL; i++) - { - inst->noiseEstLogQuantile[i] = 2048; // Q8 - inst->noiseEstDensity[i] = 153; // Q9 - } - for (i = 0; i < SIMULT; i++) - { - inst->noiseEstCounter[i] = (WebRtc_Word16)(END_STARTUP_LONG * (i + 1)) / SIMULT; - } - - // Initialize suppression filter with ones - WebRtcSpl_MemSetW16((WebRtc_Word16*)inst->noiseSupFilter, 16384, HALF_ANAL_BLOCKL); - - // Set the aggressiveness: default - inst->aggrMode = 0; - - //initialize variables for new method - inst->priorNonSpeechProb = 8192; // Q14(0.5) prior probability for speech/noise - for (i = 0; i < HALF_ANAL_BLOCKL; i++) - { - inst->prevMagnU16[i] = 0; - inst->prevNoiseU32[i] = 0; //previous noise-spectrum - inst->logLrtTimeAvgW32[i] = 0; //smooth LR ratio - inst->avgMagnPause[i] = 0; //conservative noise spectrum estimate - inst->initMagnEst[i] = 0; //initial average magnitude spectrum - } - - //feature quantities - inst->thresholdSpecDiff = 50; //threshold for difference feature: determined on-line - inst->thresholdSpecFlat = 20480; //threshold for flatness: determined on-line - inst->featureLogLrt = inst->thresholdLogLrt; //average LRT factor (= threshold) - inst->featureSpecFlat = inst->thresholdSpecFlat; //spectral flatness (= threshold) - inst->featureSpecDiff = inst->thresholdSpecDiff; //spectral difference (= threshold) - inst->weightLogLrt = 6; //default weighting par for LRT feature - inst->weightSpecFlat = 0; //default weighting par for spectral flatness feature - inst->weightSpecDiff = 0; //default weighting par for spectral difference feature - - inst->curAvgMagnEnergy = 0; //window time-average of input magnitude spectrum - inst->timeAvgMagnEnergy = 0; //normalization for spectral difference - inst->timeAvgMagnEnergyTmp = 0; //normalization for spectral difference - - //histogram quantities: used to estimate/update thresholds for features - WebRtcSpl_ZerosArrayW16(inst->histLrt, HIST_PAR_EST); - WebRtcSpl_ZerosArrayW16(inst->histSpecDiff, HIST_PAR_EST); - WebRtcSpl_ZerosArrayW16(inst->histSpecFlat, HIST_PAR_EST); - - inst->blockIndex = -1; //frame counter - - //inst->modelUpdate = 500; //window for update - inst->modelUpdate = (1 << STAT_UPDATES); //window for update - inst->cntThresUpdate = 0; //counter feature thresholds updates - - inst->sumMagn = 0; - inst->magnEnergy = 0; - inst->prevQMagn = 0; - inst->qNoise = 0; - inst->prevQNoise = 0; - - inst->energyIn = 0; - inst->scaleEnergyIn = 0; - - inst->whiteNoiseLevel = 0; - inst->pinkNoiseNumerator = 0; - inst->pinkNoiseExp = 0; - inst->minNorm = 15; // Start with full scale - inst->zeroInputSignal = 0; - - //default mode - WebRtcNsx_set_policy_core(inst, 0); - -#ifdef NS_FILEDEBUG - inst->infile=fopen("indebug.pcm","wb"); - inst->outfile=fopen("outdebug.pcm","wb"); - inst->file1=fopen("file1.pcm","wb"); - inst->file2=fopen("file2.pcm","wb"); - inst->file3=fopen("file3.pcm","wb"); - inst->file4=fopen("file4.pcm","wb"); - inst->file5=fopen("file5.pcm","wb"); -#endif - - inst->initFlag = 1; - - return 0; -} - -int WebRtcNsx_set_policy_core(NsxInst_t *inst, int mode) -{ - // allow for modes:0,1,2,3 - if (mode < 0 || mode > 3) - { - return -1; - } - - inst->aggrMode = mode; - if (mode == 0) - { - inst->overdrive = 256; // Q8(1.0) - inst->denoiseBound = 8192; // Q14(0.5) - inst->gainMap = 0; // No gain compensation - } else if (mode == 1) - { - inst->overdrive = 256; // Q8(1.0) - inst->denoiseBound = 4096; // Q14(0.25) - inst->factor2Table = kFactor2Aggressiveness1; - inst->gainMap = 1; - } else if (mode == 2) - { - inst->overdrive = 282; // ~= Q8(1.1) - inst->denoiseBound = 2048; // Q14(0.125) - inst->factor2Table = kFactor2Aggressiveness2; - inst->gainMap = 1; - } else if (mode == 3) - { - inst->overdrive = 320; // Q8(1.25) - inst->denoiseBound = 1475; // ~= Q14(0.09) - inst->factor2Table = kFactor2Aggressiveness3; - inst->gainMap = 1; - } - return 0; -} - -void WebRtcNsx_NoiseEstimation(NsxInst_t *inst, WebRtc_UWord16 *magn, WebRtc_UWord32 *noise, - WebRtc_Word16 *qNoise) -{ - WebRtc_Word32 numerator; - - WebRtc_Word16 lmagn[HALF_ANAL_BLOCKL], counter, countDiv, countProd, delta, zeros, frac; - WebRtc_Word16 log2, tabind, logval, tmp16, tmp16no1, tmp16no2; - WebRtc_Word16 log2Const = 22713; // Q15 - WebRtc_Word16 widthFactor = 21845; - - int i, s, offset; - - numerator = FACTOR_Q16; - - tabind = inst->stages - inst->normData; - if (tabind < 0) - { - logval = -kLogTable[-tabind]; - } else - { - logval = kLogTable[tabind]; - } - - // lmagn(i)=log(magn(i))=log(2)*log2(magn(i)) - // magn is in Q(-stages), and the real lmagn values are: - // real_lmagn(i)=log(magn(i)*2^stages)=log(magn(i))+log(2^stages) - // lmagn in Q8 - for (i = 0; i < inst->magnLen; i++) - { - if (magn[i]) - { - zeros = WebRtcSpl_NormU32((WebRtc_UWord32)magn[i]); - frac = (WebRtc_Word16)((((WebRtc_UWord32)magn[i] << zeros) & 0x7FFFFFFF) >> 23); - // log2(magn(i)) - log2 = (WebRtc_Word16)(((31 - zeros) << 8) + kLogTableFrac[frac]); - // log2(magn(i))*log(2) - lmagn[i] = (WebRtc_Word16)WEBRTC_SPL_MUL_16_16_RSFT(log2, log2Const, 15); - // + log(2^stages) - lmagn[i] += logval; - } else - { - lmagn[i] = logval;//0; - } - } - - // loop over simultaneous estimates - for (s = 0; s < SIMULT; s++) - { - offset = s * inst->magnLen; - - // Get counter values from state - counter = inst->noiseEstCounter[s]; - countDiv = kCounterDiv[counter]; - countProd = (WebRtc_Word16)WEBRTC_SPL_MUL_16_16(counter, countDiv); - - // quant_est(...) - for (i = 0; i < inst->magnLen; i++) - { - // compute delta - if (inst->noiseEstDensity[offset + i] > 512) - { - delta = WebRtcSpl_DivW32W16ResW16(numerator, - inst->noiseEstDensity[offset + i]); - } else - { - delta = FACTOR_Q7; - } - - // update log quantile estimate - tmp16 = (WebRtc_Word16)WEBRTC_SPL_MUL_16_16_RSFT(delta, countDiv, 14); - if (lmagn[i] > inst->noiseEstLogQuantile[offset + i]) - { - // +=QUANTILE*delta/(inst->counter[s]+1) QUANTILE=0.25, =1 in Q2 - // CounterDiv=1/inst->counter[s] in Q15 - tmp16 += 2; - tmp16no1 = WEBRTC_SPL_RSHIFT_W16(tmp16, 2); - inst->noiseEstLogQuantile[offset + i] += tmp16no1; - } else - { - tmp16 += 1; - tmp16no1 = WEBRTC_SPL_RSHIFT_W16(tmp16, 1); - // *(1-QUANTILE), in Q2 QUANTILE=0.25, 1-0.25=0.75=3 in Q2 - tmp16no2 = (WebRtc_Word16)WEBRTC_SPL_MUL_16_16_RSFT(tmp16no1, 3, 1); - inst->noiseEstLogQuantile[offset + i] -= tmp16no2; - } - - // update density estimate - if (WEBRTC_SPL_ABS_W16(lmagn[i] - inst->noiseEstLogQuantile[offset + i]) - < WIDTH_Q8) - { - tmp16no1 = (WebRtc_Word16)WEBRTC_SPL_MUL_16_16_RSFT_WITH_ROUND( - inst->noiseEstDensity[offset + i], countProd, 15); - tmp16no2 = (WebRtc_Word16)WEBRTC_SPL_MUL_16_16_RSFT_WITH_ROUND(widthFactor, - countDiv, 15); - inst->noiseEstDensity[offset + i] = tmp16no1 + tmp16no2; - } - } // end loop over magnitude spectrum - - if (counter >= END_STARTUP_LONG) - { - inst->noiseEstCounter[s] = 0; - if (inst->blockIndex >= END_STARTUP_LONG) - { - WebRtcNsx_UpdateNoiseEstimate(inst, offset); - } - } - inst->noiseEstCounter[s]++; - - } // end loop over simultaneous estimates - - // Sequentially update the noise during startup - if (inst->blockIndex < END_STARTUP_LONG) - { - WebRtcNsx_UpdateNoiseEstimate(inst, offset); - } - - for (i = 0; i < inst->magnLen; i++) - { - noise[i] = (WebRtc_UWord32)(inst->noiseEstQuantile[i]); // Q(qNoise) - } - (*qNoise) = (WebRtc_Word16)inst->qNoise; -} - -// Extract thresholds for feature parameters -// histograms are computed over some window_size (given by window_pars) -// thresholds and weights are extracted every window -// flag 0 means update histogram only, flag 1 means compute the thresholds/weights -// threshold and weights are returned in: inst->priorModelPars -void WebRtcNsx_FeatureParameterExtraction(NsxInst_t *inst, int flag) -{ - WebRtc_UWord32 tmpU32; - WebRtc_UWord32 histIndex; - WebRtc_UWord32 posPeak1SpecFlatFX, posPeak2SpecFlatFX; - WebRtc_UWord32 posPeak1SpecDiffFX, posPeak2SpecDiffFX; - - WebRtc_Word32 tmp32; - WebRtc_Word32 fluctLrtFX, thresFluctLrtFX; - WebRtc_Word32 avgHistLrtFX, avgSquareHistLrtFX, avgHistLrtComplFX; - - WebRtc_Word16 j; - WebRtc_Word16 numHistLrt; - - int i; - int useFeatureSpecFlat, useFeatureSpecDiff, featureSum; - int maxPeak1, maxPeak2; - int weightPeak1SpecFlat, weightPeak2SpecFlat; - int weightPeak1SpecDiff, weightPeak2SpecDiff; - - //update histograms - if (!flag) - { - // LRT - // Type casting to UWord32 is safe since negative values will not be wrapped to larger - // values than HIST_PAR_EST - histIndex = (WebRtc_UWord32)(inst->featureLogLrt); - if (histIndex < HIST_PAR_EST) - { - inst->histLrt[histIndex]++; - } - // Spectral flatness - // (inst->featureSpecFlat*20)>>10 = (inst->featureSpecFlat*5)>>8 - histIndex = WEBRTC_SPL_RSHIFT_U32(inst->featureSpecFlat * 5, 8); - if (histIndex < HIST_PAR_EST) - { - inst->histSpecFlat[histIndex]++; - } - // Spectral difference - histIndex = HIST_PAR_EST; - if (inst->timeAvgMagnEnergy) - { - // Guard against division by zero - // If timeAvgMagnEnergy == 0 we have no normalizing statistics and therefore can't - // update the histogram - histIndex = WEBRTC_SPL_UDIV((inst->featureSpecDiff * 5) >> inst->stages, - inst->timeAvgMagnEnergy); - } - if (histIndex < HIST_PAR_EST) - { - inst->histSpecDiff[histIndex]++; - } - } - - // extract parameters for speech/noise probability - if (flag) - { - useFeatureSpecDiff = 1; - //for LRT feature: - // compute the average over inst->featureExtractionParams.rangeAvgHistLrt - avgHistLrtFX = 0; - avgSquareHistLrtFX = 0; - numHistLrt = 0; - for (i = 0; i < BIN_SIZE_LRT; i++) - { - j = (2 * i + 1); - tmp32 = WEBRTC_SPL_MUL_16_16(inst->histLrt[i], j); - avgHistLrtFX += tmp32; - numHistLrt += inst->histLrt[i]; - avgSquareHistLrtFX += WEBRTC_SPL_MUL_32_16(tmp32, j); - } - avgHistLrtComplFX = avgHistLrtFX; - for (; i < HIST_PAR_EST; i++) - { - j = (2 * i + 1); - tmp32 = WEBRTC_SPL_MUL_16_16(inst->histLrt[i], j); - avgHistLrtComplFX += tmp32; - avgSquareHistLrtFX += WEBRTC_SPL_MUL_32_16(tmp32, j); - } - fluctLrtFX = WEBRTC_SPL_MUL(avgSquareHistLrtFX, numHistLrt); - fluctLrtFX -= WEBRTC_SPL_MUL(avgHistLrtFX, avgHistLrtComplFX); - thresFluctLrtFX = THRES_FLUCT_LRT * numHistLrt; - // get threshold for LRT feature: - tmpU32 = (FACTOR_1_LRT_DIFF * (WebRtc_UWord32)avgHistLrtFX); - if ((fluctLrtFX < thresFluctLrtFX) || (numHistLrt == 0) || (tmpU32 - > (WebRtc_UWord32)(100 * numHistLrt))) - { - inst->thresholdLogLrt = inst->maxLrt; //very low fluctuation, so likely noise - } else - { - tmp32 = (WebRtc_Word32)((tmpU32 << (9 + inst->stages)) / numHistLrt / 25); - // check if value is within min/max range - inst->thresholdLogLrt = WEBRTC_SPL_SAT(inst->maxLrt, tmp32, inst->minLrt); - } - if (fluctLrtFX < thresFluctLrtFX) - { - // Do not use difference feature if fluctuation of LRT feature is very low: - // most likely just noise state - useFeatureSpecDiff = 0; - } - - // for spectral flatness and spectral difference: compute the main peaks of histogram - maxPeak1 = 0; - maxPeak2 = 0; - posPeak1SpecFlatFX = 0; - posPeak2SpecFlatFX = 0; - weightPeak1SpecFlat = 0; - weightPeak2SpecFlat = 0; - - // peaks for flatness - for (i = 0; i < HIST_PAR_EST; i++) - { - if (inst->histSpecFlat[i] > maxPeak1) - { - // Found new "first" peak - maxPeak2 = maxPeak1; - weightPeak2SpecFlat = weightPeak1SpecFlat; - posPeak2SpecFlatFX = posPeak1SpecFlatFX; - - maxPeak1 = inst->histSpecFlat[i]; - weightPeak1SpecFlat = inst->histSpecFlat[i]; - posPeak1SpecFlatFX = (WebRtc_UWord32)(2 * i + 1); - } else if (inst->histSpecFlat[i] > maxPeak2) - { - // Found new "second" peak - maxPeak2 = inst->histSpecFlat[i]; - weightPeak2SpecFlat = inst->histSpecFlat[i]; - posPeak2SpecFlatFX = (WebRtc_UWord32)(2 * i + 1); - } - } - - // for spectral flatness feature - useFeatureSpecFlat = 1; - // merge the two peaks if they are close - if ((posPeak1SpecFlatFX - posPeak2SpecFlatFX < LIM_PEAK_SPACE_FLAT_DIFF) - && (weightPeak2SpecFlat * LIM_PEAK_WEIGHT_FLAT_DIFF > weightPeak1SpecFlat)) - { - weightPeak1SpecFlat += weightPeak2SpecFlat; - posPeak1SpecFlatFX = (posPeak1SpecFlatFX + posPeak2SpecFlatFX) >> 1; - } - //reject if weight of peaks is not large enough, or peak value too small - if (weightPeak1SpecFlat < THRES_WEIGHT_FLAT_DIFF || posPeak1SpecFlatFX - < THRES_PEAK_FLAT) - { - useFeatureSpecFlat = 0; - } else // if selected, get the threshold - { - // compute the threshold and check if value is within min/max range - inst->thresholdSpecFlat = WEBRTC_SPL_SAT(MAX_FLAT_Q10, FACTOR_2_FLAT_Q10 - * posPeak1SpecFlatFX, MIN_FLAT_Q10); //Q10 - } - // done with flatness feature - - if (useFeatureSpecDiff) - { - //compute two peaks for spectral difference - maxPeak1 = 0; - maxPeak2 = 0; - posPeak1SpecDiffFX = 0; - posPeak2SpecDiffFX = 0; - weightPeak1SpecDiff = 0; - weightPeak2SpecDiff = 0; - // peaks for spectral difference - for (i = 0; i < HIST_PAR_EST; i++) - { - if (inst->histSpecDiff[i] > maxPeak1) - { - // Found new "first" peak - maxPeak2 = maxPeak1; - weightPeak2SpecDiff = weightPeak1SpecDiff; - posPeak2SpecDiffFX = posPeak1SpecDiffFX; - - maxPeak1 = inst->histSpecDiff[i]; - weightPeak1SpecDiff = inst->histSpecDiff[i]; - posPeak1SpecDiffFX = (WebRtc_UWord32)(2 * i + 1); - } else if (inst->histSpecDiff[i] > maxPeak2) - { - // Found new "second" peak - maxPeak2 = inst->histSpecDiff[i]; - weightPeak2SpecDiff = inst->histSpecDiff[i]; - posPeak2SpecDiffFX = (WebRtc_UWord32)(2 * i + 1); - } - } - - // merge the two peaks if they are close - if ((posPeak1SpecDiffFX - posPeak2SpecDiffFX < LIM_PEAK_SPACE_FLAT_DIFF) - && (weightPeak2SpecDiff * LIM_PEAK_WEIGHT_FLAT_DIFF > weightPeak1SpecDiff)) - { - weightPeak1SpecDiff += weightPeak2SpecDiff; - posPeak1SpecDiffFX = (posPeak1SpecDiffFX + posPeak2SpecDiffFX) >> 1; - } - // get the threshold value and check if value is within min/max range - inst->thresholdSpecDiff = WEBRTC_SPL_SAT(MAX_DIFF, FACTOR_1_LRT_DIFF - * posPeak1SpecDiffFX, MIN_DIFF); //5x bigger - //reject if weight of peaks is not large enough - if (weightPeak1SpecDiff < THRES_WEIGHT_FLAT_DIFF) - { - useFeatureSpecDiff = 0; - } - // done with spectral difference feature - } - - // select the weights between the features - // inst->priorModelPars[4] is weight for LRT: always selected - featureSum = 6 / (1 + useFeatureSpecFlat + useFeatureSpecDiff); - inst->weightLogLrt = featureSum; - inst->weightSpecFlat = useFeatureSpecFlat * featureSum; - inst->weightSpecDiff = useFeatureSpecDiff * featureSum; - - // set histograms to zero for next update - WebRtcSpl_ZerosArrayW16(inst->histLrt, HIST_PAR_EST); - WebRtcSpl_ZerosArrayW16(inst->histSpecDiff, HIST_PAR_EST); - WebRtcSpl_ZerosArrayW16(inst->histSpecFlat, HIST_PAR_EST); - } // end of flag == 1 -} - - -// Compute spectral flatness on input spectrum -// magn is the magnitude spectrum -// spectral flatness is returned in inst->featureSpecFlat -void WebRtcNsx_ComputeSpectralFlatness(NsxInst_t *inst, WebRtc_UWord16 *magn) -{ - WebRtc_UWord32 tmpU32; - WebRtc_UWord32 avgSpectralFlatnessNum, avgSpectralFlatnessDen; - - WebRtc_Word32 tmp32; - WebRtc_Word32 currentSpectralFlatness, logCurSpectralFlatness; - - WebRtc_Word16 zeros, frac, intPart; - - int i; - - // for flatness - avgSpectralFlatnessNum = 0; - avgSpectralFlatnessDen = inst->sumMagn - (WebRtc_UWord32)magn[0]; // Q(normData-stages) - - // compute log of ratio of the geometric to arithmetic mean: check for log(0) case - // flatness = exp( sum(log(magn[i]))/N - log(sum(magn[i])/N) ) - // = exp( sum(log(magn[i]))/N ) * N / sum(magn[i]) - // = 2^( sum(log2(magn[i]))/N - (log2(sum(magn[i])) - log2(N)) ) [This is used] - for (i = 1; i < inst->magnLen; i++) - { - // First bin is excluded from spectrum measures. Number of bins is now a power of 2 - if (magn[i]) - { - zeros = WebRtcSpl_NormU32((WebRtc_UWord32)magn[i]); - frac = (WebRtc_Word16)(((WebRtc_UWord32)((WebRtc_UWord32)(magn[i]) << zeros) - & 0x7FFFFFFF) >> 23); - // log2(magn(i)) - tmpU32 = (WebRtc_UWord32)(((31 - zeros) << 8) + kLogTableFrac[frac]); // Q8 - avgSpectralFlatnessNum += tmpU32; // Q8 - } else - { - //if at least one frequency component is zero, treat separately - tmpU32 = WEBRTC_SPL_UMUL_32_16(inst->featureSpecFlat, SPECT_FLAT_TAVG_Q14); // Q24 - inst->featureSpecFlat -= WEBRTC_SPL_RSHIFT_U32(tmpU32, 14); // Q10 - return; - } - } - //ratio and inverse log: check for case of log(0) - zeros = WebRtcSpl_NormU32(avgSpectralFlatnessDen); - frac = (WebRtc_Word16)(((avgSpectralFlatnessDen << zeros) & 0x7FFFFFFF) >> 23); - // log2(avgSpectralFlatnessDen) - tmp32 = (WebRtc_Word32)(((31 - zeros) << 8) + kLogTableFrac[frac]); // Q8 - logCurSpectralFlatness = (WebRtc_Word32)avgSpectralFlatnessNum; - logCurSpectralFlatness += ((WebRtc_Word32)(inst->stages - 1) << (inst->stages + 7)); // Q(8+stages-1) - logCurSpectralFlatness -= (tmp32 << (inst->stages - 1)); - logCurSpectralFlatness = WEBRTC_SPL_LSHIFT_W32(logCurSpectralFlatness, 10 - inst->stages); // Q17 - tmp32 = (WebRtc_Word32)(0x00020000 | (WEBRTC_SPL_ABS_W32(logCurSpectralFlatness) - & 0x0001FFFF)); //Q17 - intPart = -(WebRtc_Word16)WEBRTC_SPL_RSHIFT_W32(logCurSpectralFlatness, 17); - intPart += 7; // Shift 7 to get the output in Q10 (from Q17 = -17+10) - if (intPart > 0) - { - currentSpectralFlatness = WEBRTC_SPL_RSHIFT_W32(tmp32, intPart); - } else - { - currentSpectralFlatness = WEBRTC_SPL_LSHIFT_W32(tmp32, -intPart); - } - - //time average update of spectral flatness feature - tmp32 = currentSpectralFlatness - (WebRtc_Word32)inst->featureSpecFlat; // Q10 - tmp32 = WEBRTC_SPL_MUL_32_16(SPECT_FLAT_TAVG_Q14, tmp32); // Q24 - inst->featureSpecFlat = (WebRtc_UWord32)((WebRtc_Word32)inst->featureSpecFlat - + WEBRTC_SPL_RSHIFT_W32(tmp32, 14)); // Q10 - // done with flatness feature -} - - -// Compute the difference measure between input spectrum and a template/learned noise spectrum -// magn_tmp is the input spectrum -// the reference/template spectrum is inst->magn_avg_pause[i] -// returns (normalized) spectral difference in inst->featureSpecDiff -void WebRtcNsx_ComputeSpectralDifference(NsxInst_t *inst, WebRtc_UWord16 *magnIn) -{ - // This is to be calculated: - // avgDiffNormMagn = var(magnIn) - cov(magnIn, magnAvgPause)^2 / var(magnAvgPause) - - WebRtc_UWord32 tmpU32no1, tmpU32no2; - WebRtc_UWord32 varMagnUFX, varPauseUFX, avgDiffNormMagnUFX; - - WebRtc_Word32 tmp32no1, tmp32no2; - WebRtc_Word32 avgPauseFX, avgMagnFX, covMagnPauseFX; - WebRtc_Word32 maxPause, minPause; - - WebRtc_Word16 tmp16no1; - - int i, norm32, nShifts; - - avgPauseFX = 0; - maxPause = 0; - minPause = inst->avgMagnPause[0]; // Q(prevQMagn) - // compute average quantities - for (i = 0; i < inst->magnLen; i++) - { - // Compute mean of magn_pause - avgPauseFX += inst->avgMagnPause[i]; // in Q(prevQMagn) - maxPause = WEBRTC_SPL_MAX(maxPause, inst->avgMagnPause[i]); - minPause = WEBRTC_SPL_MIN(minPause, inst->avgMagnPause[i]); - } - // normalize by replacing div of "inst->magnLen" with "inst->stages-1" shifts - avgPauseFX = WEBRTC_SPL_RSHIFT_W32(avgPauseFX, inst->stages - 1); - avgMagnFX = (WebRtc_Word32)WEBRTC_SPL_RSHIFT_U32(inst->sumMagn, inst->stages - 1); - // Largest possible deviation in magnPause for (co)var calculations - tmp32no1 = WEBRTC_SPL_MAX(maxPause - avgPauseFX, avgPauseFX - minPause); - // Get number of shifts to make sure we don't get wrap around in varPause - nShifts = WEBRTC_SPL_MAX(0, 10 + inst->stages - WebRtcSpl_NormW32(tmp32no1)); - - varMagnUFX = 0; - varPauseUFX = 0; - covMagnPauseFX = 0; - for (i = 0; i < inst->magnLen; i++) - { - // Compute var and cov of magn and magn_pause - tmp16no1 = (WebRtc_Word16)((WebRtc_Word32)magnIn[i] - avgMagnFX); - tmp32no2 = inst->avgMagnPause[i] - avgPauseFX; - varMagnUFX += (WebRtc_UWord32)WEBRTC_SPL_MUL_16_16(tmp16no1, tmp16no1); // Q(2*qMagn) - tmp32no1 = WEBRTC_SPL_MUL_32_16(tmp32no2, tmp16no1); // Q(prevQMagn+qMagn) - covMagnPauseFX += tmp32no1; // Q(prevQMagn+qMagn) - tmp32no1 = WEBRTC_SPL_RSHIFT_W32(tmp32no2, nShifts); // Q(prevQMagn-minPause) - varPauseUFX += (WebRtc_UWord32)WEBRTC_SPL_MUL(tmp32no1, tmp32no1); // Q(2*(prevQMagn-minPause)) - } - //update of average magnitude spectrum: Q(-2*stages) and averaging replaced by shifts - inst->curAvgMagnEnergy += WEBRTC_SPL_RSHIFT_U32(inst->magnEnergy, 2 * inst->normData - + inst->stages - 1); - - avgDiffNormMagnUFX = varMagnUFX; // Q(2*qMagn) - if ((varPauseUFX) && (covMagnPauseFX)) - { - tmpU32no1 = (WebRtc_UWord32)WEBRTC_SPL_ABS_W32(covMagnPauseFX); // Q(prevQMagn+qMagn) - norm32 = WebRtcSpl_NormU32(tmpU32no1) - 16; - if (norm32 > 0) - { - tmpU32no1 = WEBRTC_SPL_LSHIFT_U32(tmpU32no1, norm32); // Q(prevQMagn+qMagn+norm32) - } else - { - tmpU32no1 = WEBRTC_SPL_RSHIFT_U32(tmpU32no1, -norm32); // Q(prevQMagn+qMagn+norm32) - } - tmpU32no2 = WEBRTC_SPL_UMUL(tmpU32no1, tmpU32no1); // Q(2*(prevQMagn+qMagn-norm32)) - - nShifts += norm32; - nShifts <<= 1; - if (nShifts < 0) - { - varPauseUFX >>= (-nShifts); // Q(2*(qMagn+norm32+minPause)) - nShifts = 0; - } - tmpU32no1 = WEBRTC_SPL_UDIV(tmpU32no2, varPauseUFX); // Q(2*(qMagn+norm32-16+minPause)) - tmpU32no1 = WEBRTC_SPL_RSHIFT_U32(tmpU32no1, nShifts); - - avgDiffNormMagnUFX -= WEBRTC_SPL_MIN(avgDiffNormMagnUFX, tmpU32no1); // Q(2*qMagn) - } - //normalize and compute time average update of difference feature - tmpU32no1 = WEBRTC_SPL_RSHIFT_U32(avgDiffNormMagnUFX, 2 * inst->normData); - if (inst->featureSpecDiff > tmpU32no1) - { - tmpU32no2 = WEBRTC_SPL_UMUL_32_16(inst->featureSpecDiff - tmpU32no1, - SPECT_DIFF_TAVG_Q8); // Q(8-2*stages) - inst->featureSpecDiff -= WEBRTC_SPL_RSHIFT_U32(tmpU32no2, 8); // Q(-2*stages) - } else - { - tmpU32no2 = WEBRTC_SPL_UMUL_32_16(tmpU32no1 - inst->featureSpecDiff, - SPECT_DIFF_TAVG_Q8); // Q(8-2*stages) - inst->featureSpecDiff += WEBRTC_SPL_RSHIFT_U32(tmpU32no2, 8); // Q(-2*stages) - } -} - -// Compute speech/noise probability -// speech/noise probability is returned in: probSpeechFinal -//snrLocPrior is the prior SNR for each frequency (in Q11) -//snrLocPost is the post SNR for each frequency (in Q11) -void WebRtcNsx_SpeechNoiseProb(NsxInst_t *inst, WebRtc_UWord16 *nonSpeechProbFinal, - WebRtc_UWord32 *priorLocSnr, WebRtc_UWord32 *postLocSnr) -{ - WebRtc_UWord32 zeros, num, den, tmpU32no1, tmpU32no2, tmpU32no3; - - WebRtc_Word32 invLrtFX, indPriorFX, tmp32, tmp32no1, tmp32no2, besselTmpFX32; - WebRtc_Word32 frac32, logTmp; - WebRtc_Word32 logLrtTimeAvgKsumFX; - - WebRtc_Word16 indPriorFX16; - WebRtc_Word16 tmp16, tmp16no1, tmp16no2, tmpIndFX, tableIndex, frac, intPart; - - int i, normTmp, normTmp2, nShifts; - - // compute feature based on average LR factor - // this is the average over all frequencies of the smooth log LRT - logLrtTimeAvgKsumFX = 0; - for (i = 0; i < inst->magnLen; i++) - { - besselTmpFX32 = (WebRtc_Word32)postLocSnr[i]; // Q11 - normTmp = WebRtcSpl_NormU32(postLocSnr[i]); - num = WEBRTC_SPL_LSHIFT_U32(postLocSnr[i], normTmp); // Q(11+normTmp) - if (normTmp > 10) - { - den = WEBRTC_SPL_LSHIFT_U32(priorLocSnr[i], normTmp - 11); // Q(normTmp) - } else - { - den = WEBRTC_SPL_RSHIFT_U32(priorLocSnr[i], 11 - normTmp); // Q(normTmp) - } - besselTmpFX32 -= WEBRTC_SPL_UDIV(num, den); // Q11 - - // inst->logLrtTimeAvg[i] += LRT_TAVG * (besselTmp - log(snrLocPrior) - inst->logLrtTimeAvg[i]); - // Here, LRT_TAVG = 0.5 - zeros = WebRtcSpl_NormU32(priorLocSnr[i]); - frac32 = (WebRtc_Word32)(((priorLocSnr[i] << zeros) & 0x7FFFFFFF) >> 19); - tmp32 = WEBRTC_SPL_MUL(frac32, frac32); - tmp32 = WEBRTC_SPL_RSHIFT_W32(WEBRTC_SPL_MUL(tmp32, -43), 19); - tmp32 += WEBRTC_SPL_MUL_16_16_RSFT((WebRtc_Word16)frac32, 5412, 12); - frac32 = tmp32 + 37; - // tmp32 = log2(priorLocSnr[i]) - tmp32 = (WebRtc_Word32)(((31 - zeros) << 12) + frac32) - (11 << 12); // Q12 - logTmp = WEBRTC_SPL_RSHIFT_W32(WEBRTC_SPL_MUL_32_16(tmp32, 178), 8); // log2(priorLocSnr[i])*log(2) - tmp32no1 = WEBRTC_SPL_RSHIFT_W32(logTmp + inst->logLrtTimeAvgW32[i], 1); // Q12 - inst->logLrtTimeAvgW32[i] += (besselTmpFX32 - tmp32no1); // Q12 - - logLrtTimeAvgKsumFX += inst->logLrtTimeAvgW32[i]; // Q12 - } - inst->featureLogLrt = WEBRTC_SPL_RSHIFT_W32(logLrtTimeAvgKsumFX * 5, inst->stages + 10); // 5 = BIN_SIZE_LRT / 2 - // done with computation of LR factor - - // - //compute the indicator functions - // - - // average LRT feature - // FLOAT code - // indicator0 = 0.5 * (tanh(widthPrior * (logLrtTimeAvgKsum - threshPrior0)) + 1.0); - tmpIndFX = 16384; // Q14(1.0) - tmp32no1 = logLrtTimeAvgKsumFX - inst->thresholdLogLrt; // Q12 - nShifts = 7 - inst->stages; // WIDTH_PR_MAP_SHIFT - inst->stages + 5; - //use larger width in tanh map for pause regions - if (tmp32no1 < 0) - { - tmpIndFX = 0; - tmp32no1 = -tmp32no1; - //widthPrior = widthPrior * 2.0; - nShifts++; - } - tmp32no1 = WEBRTC_SPL_SHIFT_W32(tmp32no1, nShifts); // Q14 - // compute indicator function: sigmoid map - tableIndex = (WebRtc_Word16)WEBRTC_SPL_RSHIFT_W32(tmp32no1, 14); - if ((tableIndex < 16) && (tableIndex >= 0)) - { - tmp16no2 = kIndicatorTable[tableIndex]; - tmp16no1 = kIndicatorTable[tableIndex + 1] - kIndicatorTable[tableIndex]; - frac = (WebRtc_Word16)(tmp32no1 & 0x00003fff); // Q14 - tmp16no2 += (WebRtc_Word16)WEBRTC_SPL_MUL_16_16_RSFT(tmp16no1, frac, 14); - if (tmpIndFX == 0) - { - tmpIndFX = 8192 - tmp16no2; // Q14 - } else - { - tmpIndFX = 8192 + tmp16no2; // Q14 - } - } - indPriorFX = WEBRTC_SPL_MUL_16_16(inst->weightLogLrt, tmpIndFX); // 6*Q14 - - //spectral flatness feature - if (inst->weightSpecFlat) - { - tmpU32no1 = WEBRTC_SPL_UMUL(inst->featureSpecFlat, 400); // Q10 - tmpIndFX = 16384; // Q14(1.0) - //use larger width in tanh map for pause regions - tmpU32no2 = inst->thresholdSpecFlat - tmpU32no1; //Q10 - nShifts = 4; - if (inst->thresholdSpecFlat < tmpU32no1) - { - tmpIndFX = 0; - tmpU32no2 = tmpU32no1 - inst->thresholdSpecFlat; - //widthPrior = widthPrior * 2.0; - nShifts++; - } - tmp32no1 = (WebRtc_Word32)WebRtcSpl_DivU32U16(WEBRTC_SPL_LSHIFT_U32(tmpU32no2, - nShifts), 25); //Q14 - tmpU32no1 = WebRtcSpl_DivU32U16(WEBRTC_SPL_LSHIFT_U32(tmpU32no2, nShifts), 25); //Q14 - // compute indicator function: sigmoid map - // FLOAT code - // indicator1 = 0.5 * (tanh(sgnMap * widthPrior * (threshPrior1 - tmpFloat1)) + 1.0); - tableIndex = (WebRtc_Word16)WEBRTC_SPL_RSHIFT_U32(tmpU32no1, 14); - if (tableIndex < 16) - { - tmp16no2 = kIndicatorTable[tableIndex]; - tmp16no1 = kIndicatorTable[tableIndex + 1] - kIndicatorTable[tableIndex]; - frac = (WebRtc_Word16)(tmpU32no1 & 0x00003fff); // Q14 - tmp16no2 += (WebRtc_Word16)WEBRTC_SPL_MUL_16_16_RSFT(tmp16no1, frac, 14); - if (tmpIndFX) - { - tmpIndFX = 8192 + tmp16no2; // Q14 - } else - { - tmpIndFX = 8192 - tmp16no2; // Q14 - } - } - indPriorFX += WEBRTC_SPL_MUL_16_16(inst->weightSpecFlat, tmpIndFX); // 6*Q14 - } - - //for template spectral-difference - if (inst->weightSpecDiff) - { - tmpU32no1 = 0; - if (inst->featureSpecDiff) - { - normTmp = WEBRTC_SPL_MIN(20 - inst->stages, - WebRtcSpl_NormU32(inst->featureSpecDiff)); - tmpU32no1 = WEBRTC_SPL_LSHIFT_U32(inst->featureSpecDiff, normTmp); // Q(normTmp-2*stages) - tmpU32no2 = WEBRTC_SPL_RSHIFT_U32(inst->timeAvgMagnEnergy, 20 - inst->stages - - normTmp); - if (tmpU32no2) - { - tmpU32no1 = WEBRTC_SPL_UDIV(tmpU32no1, tmpU32no2); // Q14?? Q(20 - inst->stages) - } else - { - tmpU32no1 = (WebRtc_UWord32)(0x7fffffff); - } - } - tmpU32no3 = WEBRTC_SPL_UDIV(WEBRTC_SPL_LSHIFT_U32(inst->thresholdSpecDiff, 17), 25); - tmpU32no2 = tmpU32no1 - tmpU32no3; - nShifts = 1; - tmpIndFX = 16384; // Q14(1.0) - //use larger width in tanh map for pause regions - if (tmpU32no2 & 0x80000000) - { - tmpIndFX = 0; - tmpU32no2 = tmpU32no3 - tmpU32no1; - //widthPrior = widthPrior * 2.0; - nShifts--; - } - tmpU32no1 = WEBRTC_SPL_RSHIFT_U32(tmpU32no2, nShifts); - // compute indicator function: sigmoid map - /* FLOAT code - indicator2 = 0.5 * (tanh(widthPrior * (tmpFloat1 - threshPrior2)) + 1.0); - */ - tableIndex = (WebRtc_Word16)WEBRTC_SPL_RSHIFT_U32(tmpU32no1, 14); - if (tableIndex < 16) - { - tmp16no2 = kIndicatorTable[tableIndex]; - tmp16no1 = kIndicatorTable[tableIndex + 1] - kIndicatorTable[tableIndex]; - frac = (WebRtc_Word16)(tmpU32no1 & 0x00003fff); // Q14 - tmp16no2 += (WebRtc_Word16)WEBRTC_SPL_MUL_16_16_RSFT_WITH_ROUND(tmp16no1, frac, - 14); - if (tmpIndFX) - { - tmpIndFX = 8192 + tmp16no2; - } else - { - tmpIndFX = 8192 - tmp16no2; - } - } - indPriorFX += WEBRTC_SPL_MUL_16_16(inst->weightSpecDiff, tmpIndFX); // 6*Q14 - } - - //combine the indicator function with the feature weights - // FLOAT code - // indPrior = 1 - (weightIndPrior0 * indicator0 + weightIndPrior1 * indicator1 + weightIndPrior2 * indicator2); - indPriorFX16 = WebRtcSpl_DivW32W16ResW16(98307 - indPriorFX, 6); // Q14 - // done with computing indicator function - - //compute the prior probability - // FLOAT code - // inst->priorNonSpeechProb += PRIOR_UPDATE * (indPriorNonSpeech - inst->priorNonSpeechProb); - tmp16 = indPriorFX16 - inst->priorNonSpeechProb; // Q14 - inst->priorNonSpeechProb += (WebRtc_Word16)WEBRTC_SPL_MUL_16_16_RSFT(PRIOR_UPDATE_Q14, - tmp16, 14); // Q14 - - //final speech probability: combine prior model with LR factor: - for (i = 0; i < inst->magnLen; i++) - { - // FLOAT code - // invLrt = exp(inst->logLrtTimeAvg[i]); - // invLrt = inst->priorSpeechProb * invLrt; - // nonSpeechProbFinal[i] = (1.0 - inst->priorSpeechProb) / (1.0 - inst->priorSpeechProb + invLrt); - // invLrt = (1.0 - inst->priorNonSpeechProb) * invLrt; - // nonSpeechProbFinal[i] = inst->priorNonSpeechProb / (inst->priorNonSpeechProb + invLrt); - nonSpeechProbFinal[i] = 0; // Q8 - if ((inst->logLrtTimeAvgW32[i] < 65300) && (inst->priorNonSpeechProb > 0)) - { - tmp32no1 = WEBRTC_SPL_RSHIFT_W32(WEBRTC_SPL_MUL(inst->logLrtTimeAvgW32[i], 23637), - 14); // Q12 - intPart = (WebRtc_Word16)WEBRTC_SPL_RSHIFT_W32(tmp32no1, 12); - if (intPart < -8) - { - intPart = -8; - } - frac = (WebRtc_Word16)(tmp32no1 & 0x00000fff); // Q12 - // Quadratic approximation of 2^frac - tmp32no2 = WEBRTC_SPL_RSHIFT_W32(frac * frac * 44, 19); // Q12 - tmp32no2 += WEBRTC_SPL_MUL_16_16_RSFT(frac, 84, 7); // Q12 - invLrtFX = WEBRTC_SPL_LSHIFT_W32(1, 8 + intPart) - + WEBRTC_SPL_SHIFT_W32(tmp32no2, intPart - 4); // Q8 - - normTmp = WebRtcSpl_NormW32(invLrtFX); - normTmp2 = WebRtcSpl_NormW16((16384 - inst->priorNonSpeechProb)); - if (normTmp + normTmp2 < 15) - { - invLrtFX = WEBRTC_SPL_RSHIFT_W32(invLrtFX, 15 - normTmp2 - normTmp); // Q(normTmp+normTmp2-7) - tmp32no1 = WEBRTC_SPL_MUL_32_16(invLrtFX, (16384 - inst->priorNonSpeechProb)); // Q(normTmp+normTmp2+7) - invLrtFX = WEBRTC_SPL_SHIFT_W32(tmp32no1, 7 - normTmp - normTmp2); // Q14 - } else - { - tmp32no1 = WEBRTC_SPL_MUL_32_16(invLrtFX, (16384 - inst->priorNonSpeechProb)); // Q22 - invLrtFX = WEBRTC_SPL_RSHIFT_W32(tmp32no1, 8); // Q14 - } - - tmp32no1 = WEBRTC_SPL_LSHIFT_W32((WebRtc_Word32)inst->priorNonSpeechProb, 8); // Q22 - nonSpeechProbFinal[i] = (WebRtc_UWord16)WEBRTC_SPL_DIV(tmp32no1, - (WebRtc_Word32)inst->priorNonSpeechProb - + invLrtFX); // Q8 - if (7 - normTmp - normTmp2 > 0) - { - nonSpeechProbFinal[i] = 0; // Q8 - } - } - } -} - -// Transform input (speechFrame) to frequency domain magnitude (magnU16) -void WebRtcNsx_DataAnalysis(NsxInst_t *inst, short *speechFrame, WebRtc_UWord16 *magnU16) -{ - - WebRtc_UWord32 tmpU32no1, tmpU32no2; - - WebRtc_Word32 tmp_1_w32 = 0; - WebRtc_Word32 tmp_2_w32 = 0; - WebRtc_Word32 sum_log_magn = 0; - WebRtc_Word32 sum_log_i_log_magn = 0; - - WebRtc_UWord16 sum_log_magn_u16 = 0; - WebRtc_UWord16 tmp_u16 = 0; - - WebRtc_Word16 sum_log_i = 0; - WebRtc_Word16 sum_log_i_square = 0; - WebRtc_Word16 frac = 0; - WebRtc_Word16 log2 = 0; - WebRtc_Word16 matrix_determinant = 0; - WebRtc_Word16 winData[ANAL_BLOCKL_MAX], maxWinData; - WebRtc_Word16 realImag[ANAL_BLOCKL_MAX << 1]; - - int i, j; - int outCFFT; - int zeros; - int net_norm = 0; - int right_shifts_in_magnU16 = 0; - int right_shifts_in_initMagnEst = 0; - - // For lower band do all processing - // update analysis buffer for L band - WEBRTC_SPL_MEMCPY_W16(inst->analysisBuffer, inst->analysisBuffer + inst->blockLen10ms, - inst->anaLen - inst->blockLen10ms); - WEBRTC_SPL_MEMCPY_W16(inst->analysisBuffer + inst->anaLen - inst->blockLen10ms, - speechFrame, inst->blockLen10ms); - - // Window data before FFT - for (i = 0; i < inst->anaLen; i++) - { - winData[i] = (WebRtc_Word16)WEBRTC_SPL_MUL_16_16_RSFT_WITH_ROUND(inst->window[i], - inst->analysisBuffer[i], - 14); // Q0 - } - // Get input energy - inst->energyIn = WebRtcSpl_Energy(winData, (int)inst->anaLen, &(inst->scaleEnergyIn)); - - // Reset zero input flag - inst->zeroInputSignal = 0; - // Acquire norm for winData - maxWinData = WebRtcSpl_MaxAbsValueW16(winData, inst->anaLen); - inst->normData = WebRtcSpl_NormW16(maxWinData); - if (maxWinData == 0) - { - // Treat zero input separately. - inst->zeroInputSignal = 1; - return; - } - - // Determine the net normalization in the frequency domain - net_norm = inst->stages - inst->normData; - // Track lowest normalization factor and use it to prevent wrap around in shifting - right_shifts_in_magnU16 = inst->normData - inst->minNorm; - right_shifts_in_initMagnEst = WEBRTC_SPL_MAX(-right_shifts_in_magnU16, 0); - inst->minNorm -= right_shifts_in_initMagnEst; - right_shifts_in_magnU16 = WEBRTC_SPL_MAX(right_shifts_in_magnU16, 0); - - // create realImag as winData interleaved with zeros (= imag. part), normalize it - for (i = 0; i < inst->anaLen; i++) - { - j = WEBRTC_SPL_LSHIFT_W16(i, 1); - realImag[j] = WEBRTC_SPL_LSHIFT_W16(winData[i], inst->normData); // Q(normData) - realImag[j + 1] = 0; // Insert zeros in imaginary part - } - - // bit-reverse position of elements in array and FFT the array - WebRtcSpl_ComplexBitReverse(realImag, inst->stages); // Q(normData-stages) - outCFFT = WebRtcSpl_ComplexFFT(realImag, inst->stages, 1); - - inst->imag[0] = 0; // Q(normData-stages) - inst->imag[inst->anaLen2] = 0; - inst->real[0] = realImag[0]; // Q(normData-stages) - inst->real[inst->anaLen2] = realImag[inst->anaLen]; - // Q(2*(normData-stages)) - inst->magnEnergy = (WebRtc_UWord32)WEBRTC_SPL_MUL_16_16(inst->real[0], inst->real[0]); - inst->magnEnergy += (WebRtc_UWord32)WEBRTC_SPL_MUL_16_16(inst->real[inst->anaLen2], - inst->real[inst->anaLen2]); - magnU16[0] = (WebRtc_UWord16)WEBRTC_SPL_ABS_W16(inst->real[0]); // Q(normData-stages) - magnU16[inst->anaLen2] = (WebRtc_UWord16)WEBRTC_SPL_ABS_W16(inst->real[inst->anaLen2]); - inst->sumMagn = (WebRtc_UWord32)magnU16[0]; // Q(normData-stages) - inst->sumMagn += (WebRtc_UWord32)magnU16[inst->anaLen2]; - - // Gather information during startup for noise parameter estimation - if (inst->blockIndex < END_STARTUP_SHORT) - { - // Switch initMagnEst to Q(minNorm-stages) - inst->initMagnEst[0] = WEBRTC_SPL_RSHIFT_U32(inst->initMagnEst[0], - right_shifts_in_initMagnEst); - inst->initMagnEst[inst->anaLen2] = - WEBRTC_SPL_RSHIFT_U32(inst->initMagnEst[inst->anaLen2], - right_shifts_in_initMagnEst); // Q(minNorm-stages) - - // Shift magnU16 to same domain as initMagnEst - tmpU32no1 = WEBRTC_SPL_RSHIFT_W32((WebRtc_UWord32)magnU16[0], - right_shifts_in_magnU16); // Q(minNorm-stages) - tmpU32no2 = WEBRTC_SPL_RSHIFT_W32((WebRtc_UWord32)magnU16[inst->anaLen2], - right_shifts_in_magnU16); // Q(minNorm-stages) - - // Update initMagnEst - inst->initMagnEst[0] += tmpU32no1; // Q(minNorm-stages) - inst->initMagnEst[inst->anaLen2] += tmpU32no2; // Q(minNorm-stages) - - log2 = 0; - if (magnU16[inst->anaLen2]) - { - // Calculate log2(magnU16[inst->anaLen2]) - zeros = WebRtcSpl_NormU32((WebRtc_UWord32)magnU16[inst->anaLen2]); - frac = (WebRtc_Word16)((((WebRtc_UWord32)magnU16[inst->anaLen2] << zeros) & - 0x7FFFFFFF) >> 23); // Q8 - // log2(magnU16(i)) in Q8 - log2 = (WebRtc_Word16)(((31 - zeros) << 8) + kLogTableFrac[frac]); - } - - sum_log_magn = (WebRtc_Word32)log2; // Q8 - // sum_log_i_log_magn in Q17 - sum_log_i_log_magn = (WEBRTC_SPL_MUL_16_16(kLogIndex[inst->anaLen2], log2) >> 3); - } - - for (i = 1; i < inst->anaLen2; i++) - { - j = WEBRTC_SPL_LSHIFT_W16(i, 1); - inst->real[i] = realImag[j]; - inst->imag[i] = -realImag[j + 1]; - // magnitude spectrum - // energy in Q(2*(normData-stages)) - tmpU32no1 = (WebRtc_UWord32)WEBRTC_SPL_MUL_16_16(realImag[j], realImag[j]); - tmpU32no1 += (WebRtc_UWord32)WEBRTC_SPL_MUL_16_16(realImag[j + 1], realImag[j + 1]); - inst->magnEnergy += tmpU32no1; // Q(2*(normData-stages)) - - magnU16[i] = (WebRtc_UWord16)WebRtcSpl_Sqrt(tmpU32no1); // Q(normData-stages) - inst->sumMagn += (WebRtc_UWord32)magnU16[i]; // Q(normData-stages) - if (inst->blockIndex < END_STARTUP_SHORT) - { - // Switch initMagnEst to Q(minNorm-stages) - inst->initMagnEst[i] = WEBRTC_SPL_RSHIFT_U32(inst->initMagnEst[i], - right_shifts_in_initMagnEst); - - // Shift magnU16 to same domain as initMagnEst, i.e., Q(minNorm-stages) - tmpU32no1 = WEBRTC_SPL_RSHIFT_W32((WebRtc_UWord32)magnU16[i], - right_shifts_in_magnU16); - // Update initMagnEst - inst->initMagnEst[i] += tmpU32no1; // Q(minNorm-stages) - - if (i >= kStartBand) - { - // For pink noise estimation. Collect data neglecting lower frequency band - log2 = 0; - if (magnU16[i]) - { - zeros = WebRtcSpl_NormU32((WebRtc_UWord32)magnU16[i]); - frac = (WebRtc_Word16)((((WebRtc_UWord32)magnU16[i] << zeros) & - 0x7FFFFFFF) >> 23); - // log2(magnU16(i)) in Q8 - log2 = (WebRtc_Word16)(((31 - zeros) << 8) + kLogTableFrac[frac]); - } - sum_log_magn += (WebRtc_Word32)log2; // Q8 - // sum_log_i_log_magn in Q17 - sum_log_i_log_magn += (WEBRTC_SPL_MUL_16_16(kLogIndex[i], log2) >> 3); - } - } - } - - //compute simplified noise model during startup - if (inst->blockIndex < END_STARTUP_SHORT) - { - // Estimate White noise - // Switch whiteNoiseLevel to Q(minNorm-stages) - inst->whiteNoiseLevel = WEBRTC_SPL_RSHIFT_U32(inst->whiteNoiseLevel, - right_shifts_in_initMagnEst); - - // Update the average magnitude spectrum, used as noise estimate. - tmpU32no1 = WEBRTC_SPL_UMUL_32_16(inst->sumMagn, inst->overdrive); - tmpU32no1 = WEBRTC_SPL_RSHIFT_U32(tmpU32no1, inst->stages + 8); - - // Replacing division above with 'stages' shifts - // Shift to same Q-domain as whiteNoiseLevel - tmpU32no1 = WEBRTC_SPL_RSHIFT_U32(tmpU32no1, right_shifts_in_magnU16); - // This operation is safe from wrap around as long as END_STARTUP_SHORT < 128 - assert(END_STARTUP_SHORT < 128); - inst->whiteNoiseLevel += tmpU32no1; // Q(minNorm-stages) - - // Estimate Pink noise parameters - // Denominator used in both parameter estimates. - // The value is only dependent on the size of the frequency band (kStartBand) - // and to reduce computational complexity stored in a table (kDeterminantEstMatrix[]) - matrix_determinant = kDeterminantEstMatrix[kStartBand]; // Q0 - sum_log_i = kSumLogIndex[kStartBand]; // Q5 - sum_log_i_square = kSumSquareLogIndex[kStartBand]; // Q2 - if (inst->fs == 8000) - { - // Adjust values to shorter blocks in narrow band. - tmp_1_w32 = (WebRtc_Word32)matrix_determinant; - tmp_1_w32 += WEBRTC_SPL_MUL_16_16_RSFT(kSumLogIndex[65], sum_log_i, 9); - tmp_1_w32 -= WEBRTC_SPL_MUL_16_16_RSFT(kSumLogIndex[65], kSumLogIndex[65], 10); - tmp_1_w32 -= WEBRTC_SPL_LSHIFT_W32((WebRtc_Word32)sum_log_i_square, 4); - tmp_1_w32 -= WEBRTC_SPL_MUL_16_16_RSFT((WebRtc_Word16)(inst->magnLen - - kStartBand), kSumSquareLogIndex[65], 2); - matrix_determinant = (WebRtc_Word16)tmp_1_w32; - sum_log_i -= kSumLogIndex[65]; // Q5 - sum_log_i_square -= kSumSquareLogIndex[65]; // Q2 - } - - // Necessary number of shifts to fit sum_log_magn in a word16 - zeros = 16 - WebRtcSpl_NormW32(sum_log_magn); - if (zeros < 0) - { - zeros = 0; - } - tmp_1_w32 = WEBRTC_SPL_LSHIFT_W32(sum_log_magn, 1); // Q9 - sum_log_magn_u16 = (WebRtc_UWord16)WEBRTC_SPL_RSHIFT_W32(tmp_1_w32, zeros);//Q(9-zeros) - - // Calculate and update pinkNoiseNumerator. Result in Q11. - tmp_2_w32 = WEBRTC_SPL_MUL_16_U16(sum_log_i_square, sum_log_magn_u16); // Q(11-zeros) - tmpU32no1 = WEBRTC_SPL_RSHIFT_U32((WebRtc_UWord32)sum_log_i_log_magn, 12); // Q5 - - // Shift the largest value of sum_log_i and tmp32no3 before multiplication - tmp_u16 = WEBRTC_SPL_LSHIFT_U16((WebRtc_UWord16)sum_log_i, 1); // Q6 - if ((WebRtc_UWord32)sum_log_i > tmpU32no1) - { - tmp_u16 = WEBRTC_SPL_RSHIFT_U16(tmp_u16, zeros); - } - else - { - tmpU32no1 = WEBRTC_SPL_RSHIFT_U32(tmpU32no1, zeros); - } - tmp_2_w32 -= (WebRtc_Word32)WEBRTC_SPL_UMUL_32_16(tmpU32no1, tmp_u16); // Q(11-zeros) - matrix_determinant = WEBRTC_SPL_RSHIFT_W16(matrix_determinant, zeros); // Q(-zeros) - tmp_2_w32 = WebRtcSpl_DivW32W16(tmp_2_w32, matrix_determinant); // Q11 - tmp_2_w32 += WEBRTC_SPL_LSHIFT_W32((WebRtc_Word32)net_norm, 11); // Q11 - if (tmp_2_w32 < 0) - { - tmp_2_w32 = 0; - } - inst->pinkNoiseNumerator += tmp_2_w32; // Q11 - - // Calculate and update pinkNoiseExp. Result in Q14. - tmp_2_w32 = WEBRTC_SPL_MUL_16_U16(sum_log_i, sum_log_magn_u16); // Q(14-zeros) - tmp_1_w32 = WEBRTC_SPL_RSHIFT_W32(sum_log_i_log_magn, 3 + zeros); - tmp_1_w32 = WEBRTC_SPL_MUL((WebRtc_Word32)(inst->magnLen - kStartBand), - tmp_1_w32); - tmp_2_w32 -= tmp_1_w32; // Q(14-zeros) - if (tmp_2_w32 > 0) - { - // If the exponential parameter is negative force it to zero, which means a - // flat spectrum. - tmp_1_w32 = WebRtcSpl_DivW32W16(tmp_2_w32, matrix_determinant); // Q14 - inst->pinkNoiseExp += WEBRTC_SPL_SAT(16384, tmp_1_w32, 0); // Q14 - } - } -} - -void WebRtcNsx_DataSynthesis(NsxInst_t *inst, short *outFrame) -{ - WebRtc_Word32 tmp32no1; - WebRtc_Word32 energyOut; - - WebRtc_Word16 realImag[ANAL_BLOCKL_MAX << 1]; - WebRtc_Word16 tmp16no1, tmp16no2; - WebRtc_Word16 energyRatio; - WebRtc_Word16 gainFactor, gainFactor1, gainFactor2; - - int i, j; - int outCIFFT; - int scaleEnergyOut = 0; - - if (inst->zeroInputSignal) - { - // synthesize the special case of zero input - // read out fully processed segment - for (i = 0; i < inst->blockLen10ms; i++) - { - outFrame[i] = inst->synthesisBuffer[i]; // Q0 - } - // update synthesis buffer - WEBRTC_SPL_MEMCPY_W16(inst->synthesisBuffer, - inst->synthesisBuffer + inst->blockLen10ms, - inst->anaLen - inst->blockLen10ms); - WebRtcSpl_ZerosArrayW16(inst->synthesisBuffer + inst->anaLen - inst->blockLen10ms, - inst->blockLen10ms); - return; - } - // Filter the data in the frequency domain - for (i = 0; i < inst->magnLen; i++) - { - inst->real[i] = (WebRtc_Word16)WEBRTC_SPL_MUL_16_16_RSFT(inst->real[i], - (WebRtc_Word16)(inst->noiseSupFilter[i]), 14); // Q(normData-stages) - inst->imag[i] = (WebRtc_Word16)WEBRTC_SPL_MUL_16_16_RSFT(inst->imag[i], - (WebRtc_Word16)(inst->noiseSupFilter[i]), 14); // Q(normData-stages) - } - // back to time domain - // Create spectrum - realImag[0] = inst->real[0]; - realImag[1] = -inst->imag[0]; - for (i = 1; i < inst->anaLen2; i++) - { - j = WEBRTC_SPL_LSHIFT_W16(i, 1); - tmp16no1 = (inst->anaLen << 1) - j; - realImag[j] = inst->real[i]; - realImag[j + 1] = -inst->imag[i]; - realImag[tmp16no1] = inst->real[i]; - realImag[tmp16no1 + 1] = inst->imag[i]; - } - realImag[inst->anaLen] = inst->real[inst->anaLen2]; - realImag[inst->anaLen + 1] = -inst->imag[inst->anaLen2]; - - // bit-reverse position of elements in array and IFFT it - WebRtcSpl_ComplexBitReverse(realImag, inst->stages); - outCIFFT = WebRtcSpl_ComplexIFFT(realImag, inst->stages, 1); - - for (i = 0; i < inst->anaLen; i++) - { - j = WEBRTC_SPL_LSHIFT_W16(i, 1); - tmp32no1 = WEBRTC_SPL_SHIFT_W32((WebRtc_Word32)realImag[j], outCIFFT - inst->normData); - inst->real[i] = (WebRtc_Word16)WEBRTC_SPL_SAT(WEBRTC_SPL_WORD16_MAX, tmp32no1, - WEBRTC_SPL_WORD16_MIN); - } - - //scale factor: only do it after END_STARTUP_LONG time - gainFactor = 8192; // 8192 = Q13(1.0) - if (inst->gainMap == 1 && - inst->blockIndex > END_STARTUP_LONG && - inst->energyIn > 0) - { - energyOut = WebRtcSpl_Energy(inst->real, (int)inst->anaLen, &scaleEnergyOut); // Q(-scaleEnergyOut) - if (scaleEnergyOut == 0 && !(energyOut & 0x7f800000)) - { - energyOut = WEBRTC_SPL_SHIFT_W32(energyOut, 8 + scaleEnergyOut - - inst->scaleEnergyIn); - } else - { - inst->energyIn = WEBRTC_SPL_RSHIFT_W32(inst->energyIn, 8 + scaleEnergyOut - - inst->scaleEnergyIn); // Q(-8-scaleEnergyOut) - } - - assert(inst->energyIn > 0); - energyRatio = (WebRtc_Word16)WEBRTC_SPL_DIV(energyOut - + WEBRTC_SPL_RSHIFT_W32(inst->energyIn, 1), inst->energyIn); // Q8 - - // // original FLOAT code - // if (gain > blim) { - // factor1=1.0+1.3*(gain-blim); - // if (gain*factor1 > 1.0) { // FLOAT - // factor1 = 1.0/gain; // FLOAT - // } - // } - // else { - // factor1=1.0; // FLOAT - // } - // - // if (gain > blim) { - // factor2=1.0; //FLOAT - // } - // else { - // //don't reduce scale too much for pause regions: attenuation here should be controlled by flooring - // factor2=1.0-0.3*(blim-gain); // FLOAT - // if (gain <= inst->denoiseBound) { - // factor2=1.0-0.3*(blim-inst->denoiseBound); // FLOAT - // } - // } - - // all done in lookup tables now - gainFactor1 = kFactor1Table[energyRatio]; // Q8 - gainFactor2 = inst->factor2Table[energyRatio]; // Q8 - - //combine both scales with speech/noise prob: note prior (priorSpeechProb) is not frequency dependent - - // factor = inst->priorSpeechProb*factor1 + (1.0-inst->priorSpeechProb)*factor2; // original code - tmp16no1 = (WebRtc_Word16)WEBRTC_SPL_MUL_16_16_RSFT(16384 - inst->priorNonSpeechProb, - gainFactor1, 14); // Q13 16384 = Q14(1.0) - tmp16no2 = (WebRtc_Word16)WEBRTC_SPL_MUL_16_16_RSFT(inst->priorNonSpeechProb, - gainFactor2, 14); // Q13; - gainFactor = tmp16no1 + tmp16no2; // Q13 - } // out of flag_gain_map==1 - - // synthesis - for (i = 0; i < inst->anaLen; i++) - { - tmp16no1 = (WebRtc_Word16)WEBRTC_SPL_MUL_16_16_RSFT_WITH_ROUND(inst->window[i], - inst->real[i], 14); // Q0, window in Q14 - tmp32no1 = WEBRTC_SPL_MUL_16_16_RSFT_WITH_ROUND(tmp16no1, gainFactor, 13); // Q0 - // Down shift with rounding - tmp16no2 = (WebRtc_Word16)WEBRTC_SPL_SAT(WEBRTC_SPL_WORD16_MAX, tmp32no1, - WEBRTC_SPL_WORD16_MIN); // Q0 - inst->synthesisBuffer[i] = WEBRTC_SPL_ADD_SAT_W16(inst->synthesisBuffer[i], tmp16no2); // Q0 - } - - // read out fully processed segment - for (i = 0; i < inst->blockLen10ms; i++) - { - outFrame[i] = inst->synthesisBuffer[i]; // Q0 - } - // update synthesis buffer - WEBRTC_SPL_MEMCPY_W16(inst->synthesisBuffer, inst->synthesisBuffer + inst->blockLen10ms, - inst->anaLen - inst->blockLen10ms); - WebRtcSpl_ZerosArrayW16(inst->synthesisBuffer + inst->anaLen - inst->blockLen10ms, - inst->blockLen10ms); -} - -int WebRtcNsx_ProcessCore(NsxInst_t *inst, short *speechFrame, short *speechFrameHB, - short *outFrame, short *outFrameHB) -{ - // main routine for noise suppression - - WebRtc_UWord32 tmpU32no1, tmpU32no2, tmpU32no3; - WebRtc_UWord32 satMax, maxNoiseU32; - WebRtc_UWord32 tmpMagnU32, tmpNoiseU32; - WebRtc_UWord32 nearMagnEst; - WebRtc_UWord32 noiseUpdateU32; - WebRtc_UWord32 noiseU32[HALF_ANAL_BLOCKL]; - WebRtc_UWord32 postLocSnr[HALF_ANAL_BLOCKL]; - WebRtc_UWord32 priorLocSnr[HALF_ANAL_BLOCKL]; - WebRtc_UWord32 prevNearSnr[HALF_ANAL_BLOCKL]; - WebRtc_UWord32 curNearSnr; - WebRtc_UWord32 priorSnr; - WebRtc_UWord32 noise_estimate = 0; - WebRtc_UWord32 noise_estimate_avg = 0; - WebRtc_UWord32 numerator = 0; - - WebRtc_Word32 tmp32no1, tmp32no2; - WebRtc_Word32 pink_noise_num_avg = 0; - - WebRtc_UWord16 tmpU16no1; - WebRtc_UWord16 magnU16[HALF_ANAL_BLOCKL]; - WebRtc_UWord16 prevNoiseU16[HALF_ANAL_BLOCKL]; - WebRtc_UWord16 nonSpeechProbFinal[HALF_ANAL_BLOCKL]; - WebRtc_UWord16 gammaNoise, prevGammaNoise; - WebRtc_UWord16 noiseSupFilterTmp[HALF_ANAL_BLOCKL]; - - WebRtc_Word16 qMagn, qNoise; - WebRtc_Word16 avgProbSpeechHB, gainModHB, avgFilterGainHB, gainTimeDomainHB; - WebRtc_Word16 tmp16no1; - WebRtc_Word16 int_part = 0; - WebRtc_Word16 frac_part = 0; - WebRtc_Word16 pink_noise_exp_avg = 0; - - int i; - int nShifts, postShifts; - int norm32no1, norm32no2; - int flag, sign; - int q_domain_to_use = 0; - -#ifdef NS_FILEDEBUG - fwrite(spframe, sizeof(short), inst->blockLen10ms, inst->infile); -#endif - - // Check that initialization has been done - if (inst->initFlag != 1) - { - return -1; - } - // Check for valid pointers based on sampling rate - if ((inst->fs == 32000) && (speechFrameHB == NULL)) - { - return -1; - } - - // Store speechFrame and transform to frequency domain - WebRtcNsx_DataAnalysis(inst, speechFrame, magnU16); - - if (inst->zeroInputSignal) - { - WebRtcNsx_DataSynthesis(inst, outFrame); - - if (inst->fs == 32000) - { - // update analysis buffer for H band - // append new data to buffer FX - WEBRTC_SPL_MEMCPY_W16(inst->dataBufHBFX, inst->dataBufHBFX + inst->blockLen10ms, - inst->anaLen - inst->blockLen10ms); - WEBRTC_SPL_MEMCPY_W16(inst->dataBufHBFX + inst->anaLen - inst->blockLen10ms, - speechFrameHB, inst->blockLen10ms); - for (i = 0; i < inst->blockLen10ms; i++) - { - outFrameHB[i] = inst->dataBufHBFX[i]; // Q0 - } - } // end of H band gain computation - return 0; - } - - // Update block index when we have something to process - inst->blockIndex++; - // - - // Norm of magn - qMagn = inst->normData - inst->stages; - - // Compute spectral flatness on input spectrum - WebRtcNsx_ComputeSpectralFlatness(inst, magnU16); - - // quantile noise estimate - WebRtcNsx_NoiseEstimation(inst, magnU16, noiseU32, &qNoise); - - //noise estimate from previous frame - for (i = 0; i < inst->magnLen; i++) - { - prevNoiseU16[i] = (WebRtc_UWord16)WEBRTC_SPL_RSHIFT_U32(inst->prevNoiseU32[i], 11); // Q(prevQNoise) - } - - if (inst->blockIndex < END_STARTUP_SHORT) - { - // Noise Q-domain to be used later; see description at end of section. - q_domain_to_use = WEBRTC_SPL_MIN((int)qNoise, inst->minNorm - inst->stages); - - // Calculate frequency independent parts in parametric noise estimate and calculate - // the estimate for the lower frequency band (same values for all frequency bins) - if (inst->pinkNoiseExp) - { - pink_noise_exp_avg = (WebRtc_Word16)WebRtcSpl_DivW32W16(inst->pinkNoiseExp, - (WebRtc_Word16)(inst->blockIndex + 1)); // Q14 - pink_noise_num_avg = WebRtcSpl_DivW32W16(inst->pinkNoiseNumerator, - (WebRtc_Word16)(inst->blockIndex + 1)); // Q11 - WebRtcNsx_CalcParametricNoiseEstimate(inst, - pink_noise_exp_avg, - pink_noise_num_avg, - kStartBand, - &noise_estimate, - &noise_estimate_avg); - } - else - { - // Use white noise estimate if we have poor pink noise parameter estimates - noise_estimate = inst->whiteNoiseLevel; // Q(minNorm-stages) - noise_estimate_avg = noise_estimate / (inst->blockIndex + 1); // Q(minNorm-stages) - } - for (i = 0; i < inst->magnLen; i++) - { - // Estimate the background noise using the pink noise parameters if permitted - if ((inst->pinkNoiseExp) && (i >= kStartBand)) - { - // Reset noise_estimate - noise_estimate = 0; - noise_estimate_avg = 0; - // Calculate the parametric noise estimate for current frequency bin - WebRtcNsx_CalcParametricNoiseEstimate(inst, - pink_noise_exp_avg, - pink_noise_num_avg, - i, - &noise_estimate, - &noise_estimate_avg); - } - // Calculate parametric Wiener filter - noiseSupFilterTmp[i] = inst->denoiseBound; - if (inst->initMagnEst[i]) - { - // numerator = (initMagnEst - noise_estimate * overdrive) - // Result in Q(8+minNorm-stages) - tmpU32no1 = WEBRTC_SPL_UMUL_32_16(noise_estimate, inst->overdrive); - numerator = WEBRTC_SPL_LSHIFT_U32(inst->initMagnEst[i], 8); - if (numerator > tmpU32no1) - { - // Suppression filter coefficient larger than zero, so calculate. - numerator -= tmpU32no1; - - // Determine number of left shifts in numerator for best accuracy after - // division - nShifts = WebRtcSpl_NormU32(numerator); - nShifts = WEBRTC_SPL_SAT(6, nShifts, 0); - - // Shift numerator to Q(nShifts+8+minNorm-stages) - numerator = WEBRTC_SPL_LSHIFT_U32(numerator, nShifts); - - // Shift denominator to Q(nShifts-6+minNorm-stages) - tmpU32no1 = WEBRTC_SPL_RSHIFT_U32(inst->initMagnEst[i], 6 - nShifts); - tmpU32no2 = WEBRTC_SPL_UDIV(numerator, tmpU32no1); // Q14 - noiseSupFilterTmp[i] = (WebRtc_UWord16)WEBRTC_SPL_SAT(16384, tmpU32no2, - (WebRtc_UWord32)(inst->denoiseBound)); // Q14 - } - } - // Weight quantile noise 'noiseU32' with modeled noise 'noise_estimate_avg' - // 'noiseU32 is in Q(qNoise) and 'noise_estimate' in Q(minNorm-stages) - // To guarantee that we do not get wrap around when shifting to the same domain - // we use the lowest one. Furthermore, we need to save 6 bits for the weighting. - // 'noise_estimate_avg' can handle this operation by construction, but 'noiseU32' - // may not. - - // Shift 'noiseU32' to 'q_domain_to_use' - tmpU32no1 = WEBRTC_SPL_RSHIFT_U32(noiseU32[i], (int)qNoise - q_domain_to_use); - // Shift 'noise_estimate_avg' to 'q_domain_to_use' - tmpU32no2 = WEBRTC_SPL_RSHIFT_U32(noise_estimate_avg, inst->minNorm - inst->stages - - q_domain_to_use); - // Make a simple check to see if we have enough room for weighting 'tmpU32no1' - // without wrap around - nShifts = 0; - if (tmpU32no1 & 0xfc000000) { - tmpU32no1 = WEBRTC_SPL_RSHIFT_U32(tmpU32no1, 6); - tmpU32no2 = WEBRTC_SPL_RSHIFT_U32(tmpU32no2, 6); - nShifts = 6; - } - // Add them together and divide by startup length - noiseU32[i] = WebRtcSpl_DivU32U16(tmpU32no1 + tmpU32no2, END_STARTUP_SHORT); - // Shift back if necessary - noiseU32[i] = WEBRTC_SPL_LSHIFT_U32(noiseU32[i], nShifts); - } - // Update new Q-domain for 'noiseU32' - qNoise = q_domain_to_use; - } - // compute average signal during END_STARTUP_LONG time: - // used to normalize spectral difference measure - if (inst->blockIndex < END_STARTUP_LONG) - { - // substituting division with shift ending up in Q(-2*stages) - inst->timeAvgMagnEnergyTmp - += WEBRTC_SPL_RSHIFT_U32(inst->magnEnergy, - 2 * inst->normData + inst->stages - 1); - inst->timeAvgMagnEnergy = WebRtcSpl_DivU32U16(inst->timeAvgMagnEnergyTmp, - inst->blockIndex + 1); - } - - //start processing at frames == converged+1 - // STEP 1: compute prior and post SNR based on quantile noise estimates - - // compute direct decision (DD) estimate of prior SNR: needed for new method - satMax = (WebRtc_UWord32)1048575;// Largest possible value without getting overflow despite shifting 12 steps - postShifts = 6 + qMagn - qNoise; - nShifts = 5 - inst->prevQMagn + inst->prevQNoise; - for (i = 0; i < inst->magnLen; i++) - { - // FLOAT: - // post SNR - // postLocSnr[i] = 0.0; - // if (magn[i] > noise[i]) - // { - // postLocSnr[i] = magn[i] / (noise[i] + 0.0001); - // } - // // previous post SNR - // // previous estimate: based on previous frame with gain filter (smooth is previous filter) - // - // prevNearSnr[i] = inst->prevMagnU16[i] / (inst->noisePrev[i] + 0.0001) * (inst->smooth[i]); - // - // // DD estimate is sum of two terms: current estimate and previous estimate - // // directed decision update of priorSnr (or we actually store [2*priorSnr+1]) - // - // priorLocSnr[i] = DD_PR_SNR * prevNearSnr[i] + (1.0 - DD_PR_SNR) * (postLocSnr[i] - 1.0); - - // calculate post SNR: output in Q11 - postLocSnr[i] = 2048; // 1.0 in Q11 - tmpU32no1 = WEBRTC_SPL_LSHIFT_U32((WebRtc_UWord32)magnU16[i], 6); // Q(6+qMagn) - if (postShifts < 0) - { - tmpU32no2 = WEBRTC_SPL_RSHIFT_U32(noiseU32[i], -postShifts); // Q(6+qMagn) - } else - { - tmpU32no2 = WEBRTC_SPL_LSHIFT_U32(noiseU32[i], postShifts); // Q(6+qMagn) - } - if (tmpU32no1 > tmpU32no2) - { - // Current magnitude larger than noise - tmpU32no1 = WEBRTC_SPL_LSHIFT_U32(tmpU32no1, 11); // Q(17+qMagn) - if (tmpU32no2) - { - tmpU32no1 = WEBRTC_SPL_UDIV(tmpU32no1, tmpU32no2); // Q11 - postLocSnr[i] = WEBRTC_SPL_MIN(satMax, tmpU32no1); // Q11 - } else - { - postLocSnr[i] = satMax; - } - } - - // calculate prevNearSnr[i] and save for later instead of recalculating it later - nearMagnEst = WEBRTC_SPL_UMUL_16_16(inst->prevMagnU16[i], inst->noiseSupFilter[i]); // Q(prevQMagn+14) - tmpU32no1 = WEBRTC_SPL_LSHIFT_U32(nearMagnEst, 3); // Q(prevQMagn+17) - tmpU32no2 = WEBRTC_SPL_RSHIFT_U32(inst->prevNoiseU32[i], nShifts); // Q(prevQMagn+6) - - if (tmpU32no2) - { - tmpU32no1 = WEBRTC_SPL_DIV(tmpU32no1, tmpU32no2); // Q11 - tmpU32no1 = WEBRTC_SPL_MIN(satMax, tmpU32no1); // Q11 - } else - { - tmpU32no1 = satMax; // Q11 - } - prevNearSnr[i] = tmpU32no1; // Q11 - - //directed decision update of priorSnr - tmpU32no1 = WEBRTC_SPL_UMUL_32_16(prevNearSnr[i], DD_PR_SNR_Q11); // Q22 - tmpU32no2 = WEBRTC_SPL_UMUL_32_16(postLocSnr[i] - 2048, ONE_MINUS_DD_PR_SNR_Q11); // Q22 - priorSnr = tmpU32no1 + tmpU32no2 + 512; // Q22 (added 512 for rounding) - // priorLocSnr = 1 + 2*priorSnr - priorLocSnr[i] = 2048 + WEBRTC_SPL_RSHIFT_U32(priorSnr, 10); // Q11 - } // end of loop over frequencies - // done with step 1: DD computation of prior and post SNR - - // STEP 2: compute speech/noise likelihood - - //compute difference of input spectrum with learned/estimated noise spectrum - WebRtcNsx_ComputeSpectralDifference(inst, magnU16); - //compute histograms for determination of parameters (thresholds and weights for features) - //parameters are extracted once every window time (=inst->modelUpdate) - //counter update - inst->cntThresUpdate++; - flag = (int)(inst->cntThresUpdate == inst->modelUpdate); - //update histogram - WebRtcNsx_FeatureParameterExtraction(inst, flag); - //compute model parameters - if (flag) - { - inst->cntThresUpdate = 0; // Reset counter - //update every window: - // get normalization for spectral difference for next window estimate - - // Shift to Q(-2*stages) - inst->curAvgMagnEnergy = WEBRTC_SPL_RSHIFT_U32(inst->curAvgMagnEnergy, STAT_UPDATES); - - tmpU32no1 = (inst->curAvgMagnEnergy + inst->timeAvgMagnEnergy + 1) >> 1; //Q(-2*stages) - // Update featureSpecDiff - if ((tmpU32no1 != inst->timeAvgMagnEnergy) && (inst->featureSpecDiff)) - { - norm32no1 = 0; - tmpU32no3 = tmpU32no1; - while (0xFFFF0000 & tmpU32no3) - { - tmpU32no3 >>= 1; - norm32no1++; - } - tmpU32no2 = inst->featureSpecDiff; - while (0xFFFF0000 & tmpU32no2) - { - tmpU32no2 >>= 1; - norm32no1++; - } - tmpU32no3 = WEBRTC_SPL_UMUL(tmpU32no3, tmpU32no2); - tmpU32no3 = WEBRTC_SPL_UDIV(tmpU32no3, inst->timeAvgMagnEnergy); - if (WebRtcSpl_NormU32(tmpU32no3) < norm32no1) - { - inst->featureSpecDiff = 0x007FFFFF; - } else - { - inst->featureSpecDiff = WEBRTC_SPL_MIN(0x007FFFFF, - WEBRTC_SPL_LSHIFT_U32(tmpU32no3, norm32no1)); - } - } - - inst->timeAvgMagnEnergy = tmpU32no1; // Q(-2*stages) - inst->curAvgMagnEnergy = 0; - } - - //compute speech/noise probability - WebRtcNsx_SpeechNoiseProb(inst, nonSpeechProbFinal, priorLocSnr, postLocSnr); - - //time-avg parameter for noise update - gammaNoise = NOISE_UPDATE_Q8; // Q8 - - maxNoiseU32 = 0; - postShifts = inst->prevQNoise - qMagn; - nShifts = inst->prevQMagn - qMagn; - for (i = 0; i < inst->magnLen; i++) - { - // temporary noise update: use it for speech frames if update value is less than previous - // the formula has been rewritten into: - // noiseUpdate = noisePrev[i] + (1 - gammaNoise) * nonSpeechProb * (magn[i] - noisePrev[i]) - - if (postShifts < 0) - { - tmpU32no2 = WEBRTC_SPL_RSHIFT_U32(magnU16[i], -postShifts); // Q(prevQNoise) - } else - { - tmpU32no2 = WEBRTC_SPL_LSHIFT_U32(magnU16[i], postShifts); // Q(prevQNoise) - } - if (prevNoiseU16[i] > tmpU32no2) - { - sign = -1; - tmpU32no1 = prevNoiseU16[i] - tmpU32no2; - } else - { - sign = 1; - tmpU32no1 = tmpU32no2 - prevNoiseU16[i]; - } - noiseUpdateU32 = inst->prevNoiseU32[i]; // Q(prevQNoise+11) - tmpU32no3 = 0; - if ((tmpU32no1) && (nonSpeechProbFinal[i])) - { - // This value will be used later, if gammaNoise changes - tmpU32no3 = WEBRTC_SPL_UMUL_32_16(tmpU32no1, nonSpeechProbFinal[i]); // Q(prevQNoise+8) - if (0x7c000000 & tmpU32no3) - { - // Shifting required before multiplication - tmpU32no2 - = WEBRTC_SPL_UMUL_32_16(WEBRTC_SPL_RSHIFT_U32(tmpU32no3, 5), gammaNoise); // Q(prevQNoise+11) - } else - { - // We can do shifting after multiplication - tmpU32no2 - = WEBRTC_SPL_RSHIFT_U32(WEBRTC_SPL_UMUL_32_16(tmpU32no3, gammaNoise), 5); // Q(prevQNoise+11) - } - if (sign > 0) - { - noiseUpdateU32 += tmpU32no2; // Q(prevQNoise+11) - } else - { - // This operation is safe. We can never get wrap around, since worst - // case scenario means magnU16 = 0 - noiseUpdateU32 -= tmpU32no2; // Q(prevQNoise+11) - } - } - - //increase gamma (i.e., less noise update) for frame likely to be speech - prevGammaNoise = gammaNoise; - gammaNoise = NOISE_UPDATE_Q8; - //time-constant based on speech/noise state - //increase gamma (i.e., less noise update) for frames likely to be speech - if (nonSpeechProbFinal[i] < ONE_MINUS_PROB_RANGE_Q8) - { - gammaNoise = GAMMA_NOISE_TRANS_AND_SPEECH_Q8; - } - - if (prevGammaNoise != gammaNoise) - { - // new noise update - // this line is the same as above, only that the result is stored in a different variable and the gammaNoise - // has changed - // - // noiseUpdate = noisePrev[i] + (1 - gammaNoise) * nonSpeechProb * (magn[i] - noisePrev[i]) - - if (0x7c000000 & tmpU32no3) - { - // Shifting required before multiplication - tmpU32no2 - = WEBRTC_SPL_UMUL_32_16(WEBRTC_SPL_RSHIFT_U32(tmpU32no3, 5), gammaNoise); // Q(prevQNoise+11) - } else - { - // We can do shifting after multiplication - tmpU32no2 - = WEBRTC_SPL_RSHIFT_U32(WEBRTC_SPL_UMUL_32_16(tmpU32no3, gammaNoise), 5); // Q(prevQNoise+11) - } - if (sign > 0) - { - tmpU32no1 = inst->prevNoiseU32[i] + tmpU32no2; // Q(prevQNoise+11) - } else - { - tmpU32no1 = inst->prevNoiseU32[i] - tmpU32no2; // Q(prevQNoise+11) - } - if (noiseUpdateU32 > tmpU32no1) - { - noiseUpdateU32 = tmpU32no1; // Q(prevQNoise+11) - } - } - noiseU32[i] = noiseUpdateU32; // Q(prevQNoise+11) - if (noiseUpdateU32 > maxNoiseU32) - { - maxNoiseU32 = noiseUpdateU32; - } - - // conservative noise update - // // original FLOAT code - // if (prob_speech < PROB_RANGE) { - // inst->avgMagnPause[i] = inst->avgMagnPause[i] + (1.0 - gamma_pause)*(magn[i] - inst->avgMagnPause[i]); - // } - - tmp32no2 = WEBRTC_SPL_SHIFT_W32(inst->avgMagnPause[i], -nShifts); - if (nonSpeechProbFinal[i] > ONE_MINUS_PROB_RANGE_Q8) - { - if (nShifts < 0) - { - tmp32no1 = (WebRtc_Word32)magnU16[i] - tmp32no2; // Q(qMagn) - tmp32no1 = WEBRTC_SPL_MUL_32_16(tmp32no1, ONE_MINUS_GAMMA_PAUSE_Q8); // Q(8+prevQMagn+nShifts) - tmp32no1 = WEBRTC_SPL_RSHIFT_W32(tmp32no1 + 128, 8); // Q(qMagn) - } else - { - tmp32no1 = WEBRTC_SPL_LSHIFT_W32((WebRtc_Word32)magnU16[i], nShifts) - - inst->avgMagnPause[i]; // Q(qMagn+nShifts) - tmp32no1 = WEBRTC_SPL_MUL_32_16(tmp32no1, ONE_MINUS_GAMMA_PAUSE_Q8); // Q(8+prevQMagn+nShifts) - tmp32no1 = WEBRTC_SPL_RSHIFT_W32(tmp32no1 + (128 << nShifts), 8 + nShifts); // Q(qMagn) - } - tmp32no2 += tmp32no1; // Q(qMagn) - } - inst->avgMagnPause[i] = tmp32no2; - } // end of frequency loop - - norm32no1 = WebRtcSpl_NormU32(maxNoiseU32); - qNoise = inst->prevQNoise + norm32no1 - 5; - // done with step 2: noise update - - // STEP 3: compute dd update of prior snr and post snr based on new noise estimate - nShifts = inst->prevQNoise + 11 - qMagn; - for (i = 0; i < inst->magnLen; i++) - { - // FLOAT code - // // post and prior SNR - // curNearSnr = 0.0; - // if (magn[i] > noise[i]) - // { - // curNearSnr = magn[i] / (noise[i] + 0.0001) - 1.0; - // } - // // DD estimate is sum of two terms: current estimate and previous estimate - // // directed decision update of snrPrior - // snrPrior = DD_PR_SNR * prevNearSnr[i] + (1.0 - DD_PR_SNR) * curNearSnr; - // // gain filter - // tmpFloat1 = inst->overdrive + snrPrior; - // tmpFloat2 = snrPrior / tmpFloat1; - // theFilter[i] = tmpFloat2; - - // calculate curNearSnr again, this is necessary because a new noise estimate has been made since then. for the original - curNearSnr = 0; // Q11 - if (nShifts < 0) - { - // This case is equivalent with magn < noise which implies curNearSnr = 0; - tmpMagnU32 = (WebRtc_UWord32)magnU16[i]; // Q(qMagn) - tmpNoiseU32 = WEBRTC_SPL_LSHIFT_U32(noiseU32[i], -nShifts); // Q(qMagn) - } else if (nShifts > 17) - { - tmpMagnU32 = WEBRTC_SPL_LSHIFT_U32(magnU16[i], 17); // Q(qMagn+17) - tmpNoiseU32 = WEBRTC_SPL_RSHIFT_U32(noiseU32[i], nShifts - 17); // Q(qMagn+17) - } else - { - tmpMagnU32 = WEBRTC_SPL_LSHIFT_U32((WebRtc_UWord32)magnU16[i], nShifts); // Q(qNoise_prev+11) - tmpNoiseU32 = noiseU32[i]; // Q(qNoise_prev+11) - } - if (tmpMagnU32 > tmpNoiseU32) - { - tmpU32no1 = tmpMagnU32 - tmpNoiseU32; // Q(qCur) - norm32no2 = WEBRTC_SPL_MIN(11, WebRtcSpl_NormU32(tmpU32no1)); - tmpU32no1 = WEBRTC_SPL_LSHIFT_U32(tmpU32no1, norm32no2); // Q(qCur+norm32no2) - tmpU32no2 = WEBRTC_SPL_RSHIFT_U32(tmpNoiseU32, 11 - norm32no2); // Q(qCur+norm32no2-11) - if (tmpU32no2) - { - tmpU32no1 = WEBRTC_SPL_UDIV(tmpU32no1, tmpU32no2); // Q11 - } - curNearSnr = WEBRTC_SPL_MIN(satMax, tmpU32no1); // Q11 - } - - //directed decision update of priorSnr - // FLOAT - // priorSnr = DD_PR_SNR * prevNearSnr + (1.0-DD_PR_SNR) * curNearSnr; - - tmpU32no1 = WEBRTC_SPL_UMUL_32_16(prevNearSnr[i], DD_PR_SNR_Q11); // Q22 - tmpU32no2 = WEBRTC_SPL_UMUL_32_16(curNearSnr, ONE_MINUS_DD_PR_SNR_Q11); // Q22 - priorSnr = tmpU32no1 + tmpU32no2; // Q22 - - //gain filter - tmpU32no1 = (WebRtc_UWord32)(inst->overdrive) - + WEBRTC_SPL_RSHIFT_U32(priorSnr + 8192, 14); // Q8 - tmpU16no1 = (WebRtc_UWord16)WEBRTC_SPL_UDIV(priorSnr + (tmpU32no1 >> 1), tmpU32no1); // Q14 - inst->noiseSupFilter[i] = WEBRTC_SPL_SAT(16384, tmpU16no1, inst->denoiseBound); // 16384 = Q14(1.0) // Q14 - - // Weight in the parametric Wiener filter during startup - if (inst->blockIndex < END_STARTUP_SHORT) - { - // Weight the two suppression filters - tmpU32no1 = WEBRTC_SPL_UMUL_16_16(inst->noiseSupFilter[i], - (WebRtc_UWord16)inst->blockIndex); - tmpU32no2 = WEBRTC_SPL_UMUL_16_16(noiseSupFilterTmp[i], - (WebRtc_UWord16)(END_STARTUP_SHORT - - inst->blockIndex)); - tmpU32no1 += tmpU32no2; - inst->noiseSupFilter[i] = (WebRtc_UWord16)WebRtcSpl_DivU32U16(tmpU32no1, - END_STARTUP_SHORT); - } - } // end of loop over frequencies - //done with step3 - - // save noise and magnitude spectrum for next frame - inst->prevQNoise = qNoise; - inst->prevQMagn = qMagn; - if (norm32no1 > 5) - { - for (i = 0; i < inst->magnLen; i++) - { - inst->prevNoiseU32[i] = WEBRTC_SPL_LSHIFT_U32(noiseU32[i], norm32no1 - 5); // Q(qNoise+11) - inst->prevMagnU16[i] = magnU16[i]; // Q(qMagn) - } - } else - { - for (i = 0; i < inst->magnLen; i++) - { - inst->prevNoiseU32[i] = WEBRTC_SPL_RSHIFT_U32(noiseU32[i], 5 - norm32no1); // Q(qNoise+11) - inst->prevMagnU16[i] = magnU16[i]; // Q(qMagn) - } - } - - WebRtcNsx_DataSynthesis(inst, outFrame); -#ifdef NS_FILEDEBUG - fwrite(outframe, sizeof(short), inst->blockLen10ms, inst->outfile); -#endif - - //for H band: - // only update data buffer, then apply time-domain gain is applied derived from L band - if (inst->fs == 32000) - { - // update analysis buffer for H band - // append new data to buffer FX - WEBRTC_SPL_MEMCPY_W16(inst->dataBufHBFX, inst->dataBufHBFX + inst->blockLen10ms, inst->anaLen - inst->blockLen10ms); - WEBRTC_SPL_MEMCPY_W16(inst->dataBufHBFX + inst->anaLen - inst->blockLen10ms, speechFrameHB, inst->blockLen10ms); - // range for averaging low band quantities for H band gain - - gainTimeDomainHB = 16384; // 16384 = Q14(1.0) - //average speech prob from low band - //average filter gain from low band - //avg over second half (i.e., 4->8kHz) of freq. spectrum - tmpU32no1 = 0; // Q12 - tmpU16no1 = 0; // Q8 - for (i = inst->anaLen2 - (inst->anaLen2 >> 2); i < inst->anaLen2; i++) - { - tmpU16no1 += nonSpeechProbFinal[i]; // Q8 - tmpU32no1 += (WebRtc_UWord32)(inst->noiseSupFilter[i]); // Q14 - } - avgProbSpeechHB = (WebRtc_Word16)(4096 - - WEBRTC_SPL_RSHIFT_U16(tmpU16no1, inst->stages - 7)); // Q12 - avgFilterGainHB = (WebRtc_Word16)WEBRTC_SPL_RSHIFT_U32(tmpU32no1, inst->stages - 3); // Q14 - - // // original FLOAT code - // // gain based on speech probability: - // avg_prob_speech_tt=(float)2.0*avg_prob_speech-(float)1.0; - // gain_mod=(float)0.5*((float)1.0+(float)tanh(avg_prob_speech_tt)); // between 0 and 1 - - // gain based on speech probability: - // original expression: "0.5 * (1 + tanh(2x-1))" - // avgProbSpeechHB has been anyway saturated to a value between 0 and 1 so the other cases don't have to be dealt with - // avgProbSpeechHB and gainModHB are in Q12, 3607 = Q12(0.880615234375) which is a zero point of - // |0.5 * (1 + tanh(2x-1)) - x| - |0.5 * (1 + tanh(2x-1)) - 0.880615234375| meaning that from that point the error of approximating - // the expression with f(x) = x would be greater than the error of approximating the expression with f(x) = 0.880615234375 - // error: "|0.5 * (1 + tanh(2x-1)) - x| from x=0 to 0.880615234375" -> http://www.wolframalpha.com/input/?i=|0.5+*+(1+%2B+tanh(2x-1))+-+x|+from+x%3D0+to+0.880615234375 - // and: "|0.5 * (1 + tanh(2x-1)) - 0.880615234375| from x=0.880615234375 to 1" -> http://www.wolframalpha.com/input/?i=+|0.5+*+(1+%2B+tanh(2x-1))+-+0.880615234375|+from+x%3D0.880615234375+to+1 - gainModHB = WEBRTC_SPL_MIN(avgProbSpeechHB, 3607); - - // // original FLOAT code - // //combine gain with low band gain - // if (avg_prob_speech < (float)0.5) { - // gain_time_domain_HB=(float)0.5*gain_mod+(float)0.5*avg_filter_gain; - // } - // else { - // gain_time_domain_HB=(float)0.25*gain_mod+(float)0.75*avg_filter_gain; - // } - - - //combine gain with low band gain - if (avgProbSpeechHB < 2048) - { // 2048 = Q12(0.5) - // the next two lines in float are "gain_time_domain = 0.5 * gain_mod + 0.5 * avg_filter_gain"; Q2(0.5) = 2 equals one left shift - gainTimeDomainHB = (gainModHB << 1) + (avgFilterGainHB >> 1); // Q14 - } else - { - // "gain_time_domain = 0.25 * gain_mod + 0.75 * agv_filter_gain;" - gainTimeDomainHB = (WebRtc_Word16)WEBRTC_SPL_MUL_16_16_RSFT(3, avgFilterGainHB, 2); // 3 = Q2(0.75); Q14 - gainTimeDomainHB += gainModHB; // Q14 - } - //make sure gain is within flooring range - gainTimeDomainHB - = WEBRTC_SPL_SAT(16384, gainTimeDomainHB, (WebRtc_Word16)(inst->denoiseBound)); // 16384 = Q14(1.0) - - - //apply gain - for (i = 0; i < inst->blockLen10ms; i++) - { - outFrameHB[i] - = (WebRtc_Word16)WEBRTC_SPL_MUL_16_16_RSFT(gainTimeDomainHB, inst->dataBufHBFX[i], 14); // Q0 - } - } // end of H band gain computation - - return 0; -} |