Refactoring (bitexact) of the EchoSuppressor in WebRTC AEC (#1)

webrtc/modules/audio_processing/aec/aec_core_sse2.c

 /*
  *  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.
  */
 
@@ skipping 421 matching lines @@
   }
 }
 
 __inline static void _mm_add_ps_4x1(__m128 sum, float *dst) {
   // A+B C+D
   sum = _mm_add_ps(sum, _mm_shuffle_ps(sum, sum, _MM_SHUFFLE(0, 0, 3, 2)));
   // A+B+C+D A+B+C+D
   sum = _mm_add_ps(sum, _mm_shuffle_ps(sum, sum, _MM_SHUFFLE(1, 1, 1, 1)));
   _mm_store_ss(dst, sum);
 }
-static int PartitionDelay(const AecCore* aec) {
+
+static int PartitionDelaySSE2(const AecCore* aec) {
   // Measures the energy in each filter partition and returns the partition with
   // highest energy.
   // TODO(bjornv): Spread computational cost by computing one partition per
   // block?
   float wfEnMax = 0;
   int i;
   int delay = 0;
 
   for (i = 0; i < aec->num_partitions; i++) {
     int j;
@@ skipping 159 matching lines @@
 
   if (aec->divergeState)
     memcpy(efw, dfw, sizeof(efw[0][0]) * 2 * PART_LEN1);
 
   // Reset if error is significantly larger than nearend (13 dB).
   if (!aec->extended_filter_enabled && seSum > (19.95f * sdSum))
     memset(aec->wfBuf, 0, sizeof(aec->wfBuf));
 }
 
 // Window time domain data to be used by the fft.
-__inline static void WindowData(float* x_windowed, const float* x) {
+static void WindowDataSSE2(float* x_windowed, const float* x) {
   int i;
   for (i = 0; i < PART_LEN; i += 4) {
     const __m128 vec_Buf1 = _mm_loadu_ps(&x[i]);
     const __m128 vec_Buf2 = _mm_loadu_ps(&x[PART_LEN + i]);
     const __m128 vec_sqrtHanning = _mm_load_ps(&WebRtcAec_sqrtHanning[i]);
     // A B C D
     __m128 vec_sqrtHanning_rev =
         _mm_loadu_ps(&WebRtcAec_sqrtHanning[PART_LEN - i - 3]);
     // D C B A
     vec_sqrtHanning_rev =
         _mm_shuffle_ps(vec_sqrtHanning_rev, vec_sqrtHanning_rev,
                        _MM_SHUFFLE(0, 1, 2, 3));
     _mm_storeu_ps(&x_windowed[i], _mm_mul_ps(vec_Buf1, vec_sqrtHanning));
     _mm_storeu_ps(&x_windowed[PART_LEN + i],
                   _mm_mul_ps(vec_Buf2, vec_sqrtHanning_rev));
   }
 }
 
 // Puts fft output data into a complex valued array.
-__inline static void StoreAsComplex(const float* data,
-                                    float data_complex[2][PART_LEN1]) {
+static void StoreAsComplexSSE2(const float* data,
+                               float data_complex[2][PART_LEN1]) {
   int i;
   for (i = 0; i < PART_LEN; i += 4) {
     const __m128 vec_fft0 = _mm_loadu_ps(&data[2 * i]);
     const __m128 vec_fft4 = _mm_loadu_ps(&data[2 * i + 4]);
     const __m128 vec_a = _mm_shuffle_ps(vec_fft0, vec_fft4,
                                         _MM_SHUFFLE(2, 0, 2, 0));
     const __m128 vec_b = _mm_shuffle_ps(vec_fft0, vec_fft4,
                                         _MM_SHUFFLE(3, 1, 3, 1));
     _mm_storeu_ps(&data_complex[0][i], vec_a);
     _mm_storeu_ps(&data_complex[1][i], vec_b);
   }
   // fix beginning/end values
   data_complex[1][0] = 0;
   data_complex[1][PART_LEN] = 0;
   data_complex[0][0] = data[0];
   data_complex[0][PART_LEN] = data[1];
 }
 
 static void SubbandCoherenceSSE2(AecCore* aec,
                                  float efw[2][PART_LEN1],
+                                 float dfw[2][PART_LEN1],
                                  float xfw[2][PART_LEN1],
                                  float* fft,
                                  float* cohde,
                                  float* cohxd) {
-  float dfw[2][PART_LEN1];
   int i;
 
-  if (aec->delayEstCtr == 0)
-    aec->delayIdx = PartitionDelay(aec);
-
-  // Use delayed far.
-  memcpy(xfw,
-         aec->xfwBuf + aec->delayIdx * PART_LEN1,
-         sizeof(xfw[0][0]) * 2 * PART_LEN1);
-
-  // Windowed near fft
-  WindowData(fft, aec->dBuf);
-  aec_rdft_forward_128(fft);
-  StoreAsComplex(fft, dfw);
-
-  // Windowed error fft
-  WindowData(fft, aec->eBuf);
-  aec_rdft_forward_128(fft);
-  StoreAsComplex(fft, efw);
-
   SmoothedPSD(aec, efw, dfw, xfw);
 
   {
     const __m128 vec_1eminus10 =  _mm_set1_ps(1e-10f);
 
     // Subband coherence
     for (i = 0; i + 3 < PART_LEN1; i += 4) {
       const __m128 vec_sd = _mm_loadu_ps(&aec->sd[i]);
       const __m128 vec_se = _mm_loadu_ps(&aec->se[i]);
       const __m128 vec_sx = _mm_loadu_ps(&aec->sx[i]);
@@ skipping 34 matching lines @@
     }
   }
 }
 
 void WebRtcAec_InitAec_SSE2(void) {
   WebRtcAec_FilterFar = FilterFarSSE2;
   WebRtcAec_ScaleErrorSignal = ScaleErrorSignalSSE2;
   WebRtcAec_FilterAdaptation = FilterAdaptationSSE2;
   WebRtcAec_OverdriveAndSuppress = OverdriveAndSuppressSSE2;
   WebRtcAec_SubbandCoherence = SubbandCoherenceSSE2;
+  WebRtcAec_StoreAsComplex = StoreAsComplexSSE2;
+  WebRtcAec_PartitionDelay = PartitionDelaySSE2;
+  WebRtcAec_WindowData = WindowDataSSE2;
 }