• Home
  • History
  • Annotate
  • Line#
  • Scopes#
  • Navigate#
  • Raw
  • Download
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 
WebRtcSpl_AutoCorrelation(const int16_t * in_vector,int in_vector_length,int order,int32_t * result,int * scale)13 int WebRtcSpl_AutoCorrelation(const int16_t* in_vector,
14                               int in_vector_length,
15                               int order,
16                               int32_t* result,
17                               int* scale) {
18   int32_t sum = 0;
19   int i = 0, j = 0;
20   int16_t smax = 0;
21   int scaling = 0;
22 
23   if (order > in_vector_length) {
24     /* Undefined */
25     return -1;
26   } else if (order < 0) {
27     order = in_vector_length;
28   }
29 
30   // Find the maximum absolute value of the samples.
31   smax = WebRtcSpl_MaxAbsValueW16(in_vector, in_vector_length);
32 
33   // In order to avoid overflow when computing the sum we should scale the
34   // samples so that (in_vector_length * smax * smax) will not overflow.
35   if (smax == 0) {
36     scaling = 0;
37   } else {
38     // Number of bits in the sum loop.
39     int nbits = WebRtcSpl_GetSizeInBits(in_vector_length);
40     // Number of bits to normalize smax.
41     int t = WebRtcSpl_NormW32(WEBRTC_SPL_MUL(smax, smax));
42 
43     if (t > nbits) {
44       scaling = 0;
45     } else {
46       scaling = nbits - t;
47     }
48   }
49 
50   // Perform the actual correlation calculation.
51   for (i = 0; i < order + 1; i++) {
52     sum = 0;
53     /* Unroll the loop to improve performance. */
54     for (j = 0; j < in_vector_length - i - 3; j += 4) {
55       sum += (in_vector[j + 0] * in_vector[i + j + 0]) >> scaling;
56       sum += (in_vector[j + 1] * in_vector[i + j + 1]) >> scaling;
57       sum += (in_vector[j + 2] * in_vector[i + j + 2]) >> scaling;
58       sum += (in_vector[j + 3] * in_vector[i + j + 3]) >> scaling;
59     }
60     for (; j < in_vector_length - i; j++) {
61       sum += (in_vector[j] * in_vector[i + j]) >> scaling;
62     }
63     *result++ = sum;
64   }
65 
66   *scale = scaling;
67   return order + 1;
68 }
69