Home | History | Annotate | Download | only in signal_processing
      1 /*
      2  *  Copyright (c) 2012 The WebRTC project authors. All Rights Reserved.
      3  *
      4  *  Use of this source code is governed by a BSD-style license
      5  *  that can be found in the LICENSE file in the root of the source
      6  *  tree. An additional intellectual property rights grant can be found
      7  *  in the file PATENTS.  All contributing project authors may
      8  *  be found in the AUTHORS file in the root of the source tree.
      9  */
     10 
     11 #include "webrtc/common_audio/signal_processing/include/signal_processing_library.h"
     12 
     13 #include <assert.h>
     14 
     15 size_t WebRtcSpl_AutoCorrelation(const int16_t* in_vector,
     16                                  size_t in_vector_length,
     17                                  size_t order,
     18                                  int32_t* result,
     19                                  int* scale) {
     20   int32_t sum = 0;
     21   size_t i = 0, j = 0;
     22   int16_t smax = 0;
     23   int scaling = 0;
     24 
     25   assert(order <= in_vector_length);
     26 
     27   // Find the maximum absolute value of the samples.
     28   smax = WebRtcSpl_MaxAbsValueW16(in_vector, in_vector_length);
     29 
     30   // In order to avoid overflow when computing the sum we should scale the
     31   // samples so that (in_vector_length * smax * smax) will not overflow.
     32   if (smax == 0) {
     33     scaling = 0;
     34   } else {
     35     // Number of bits in the sum loop.
     36     int nbits = WebRtcSpl_GetSizeInBits((uint32_t)in_vector_length);
     37     // Number of bits to normalize smax.
     38     int t = WebRtcSpl_NormW32(WEBRTC_SPL_MUL(smax, smax));
     39 
     40     if (t > nbits) {
     41       scaling = 0;
     42     } else {
     43       scaling = nbits - t;
     44     }
     45   }
     46 
     47   // Perform the actual correlation calculation.
     48   for (i = 0; i < order + 1; i++) {
     49     sum = 0;
     50     /* Unroll the loop to improve performance. */
     51     for (j = 0; i + j + 3 < in_vector_length; j += 4) {
     52       sum += (in_vector[j + 0] * in_vector[i + j + 0]) >> scaling;
     53       sum += (in_vector[j + 1] * in_vector[i + j + 1]) >> scaling;
     54       sum += (in_vector[j + 2] * in_vector[i + j + 2]) >> scaling;
     55       sum += (in_vector[j + 3] * in_vector[i + j + 3]) >> scaling;
     56     }
     57     for (; j < in_vector_length - i; j++) {
     58       sum += (in_vector[j] * in_vector[i + j]) >> scaling;
     59     }
     60     *result++ = sum;
     61   }
     62 
     63   *scale = scaling;
     64   return order + 1;
     65 }
     66