Home | History | Annotate | Download | only in signal_processing 
      1 /*
      2  *  Copyright (c) 2011 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 
     12 /*
     13  * This file contains the function WebRtcSpl_AutoCorrelation().
     14  * The description header can be found in signal_processing_library.h
     15  *
     16  */
     17 
     18 #include "signal_processing_library.h"
     19 
     20 int WebRtcSpl_AutoCorrelation(G_CONST WebRtc_Word16* in_vector,
     21                               int in_vector_length,
     22                               int order,
     23                               WebRtc_Word32* result,
     24                               int* scale)
     25 {
     26     WebRtc_Word32 sum;
     27     int i, j;
     28     WebRtc_Word16 smax; // Sample max
     29     G_CONST WebRtc_Word16* xptr1;
     30     G_CONST WebRtc_Word16* xptr2;
     31     WebRtc_Word32* resultptr;
     32     int scaling = 0;
     33 
     34 #ifdef _ARM_OPT_
     35 #pragma message("NOTE: _ARM_OPT_ optimizations are used")
     36     WebRtc_Word16 loops4;
     37 #endif
     38 
     39     if (order < 0)
     40         order = in_vector_length;
     41 
     42     // Find the max. sample
     43     smax = WebRtcSpl_MaxAbsValueW16(in_vector, in_vector_length);
     44 
     45     // In order to avoid overflow when computing the sum we should scale the samples so that
     46     // (in_vector_length * smax * smax) will not overflow.
     47 
     48     if (smax == 0)
     49     {
     50         scaling = 0;
     51     } else
     52     {
     53         int nbits = WebRtcSpl_GetSizeInBits(in_vector_length); // # of bits in the sum loop
     54         int t = WebRtcSpl_NormW32(WEBRTC_SPL_MUL(smax, smax)); // # of bits to normalize smax
     55 
     56         if (t > nbits)
     57         {
     58             scaling = 0;
     59         } else
     60         {
     61             scaling = nbits - t;
     62         }
     63 
     64     }
     65 
     66     resultptr = result;
     67 
     68     // Perform the actual correlation calculation
     69     for (i = 0; i < order + 1; i++)
     70     {
     71         int loops = (in_vector_length - i);
     72         sum = 0;
     73         xptr1 = in_vector;
     74         xptr2 = &in_vector[i];
     75 #ifndef _ARM_OPT_
     76         for (j = loops; j > 0; j--)
     77         {
     78             sum += WEBRTC_SPL_MUL_16_16_RSFT(*xptr1++, *xptr2++, scaling);
     79         }
     80 #else
     81         loops4 = (loops >> 2) << 2;
     82 
     83         if (scaling == 0)
     84         {
     85             for (j = 0; j < loops4; j = j + 4)
     86             {
     87                 sum += WEBRTC_SPL_MUL_16_16(*xptr1, *xptr2);
     88                 xptr1++;
     89                 xptr2++;
     90                 sum += WEBRTC_SPL_MUL_16_16(*xptr1, *xptr2);
     91                 xptr1++;
     92                 xptr2++;
     93                 sum += WEBRTC_SPL_MUL_16_16(*xptr1, *xptr2);
     94                 xptr1++;
     95                 xptr2++;
     96                 sum += WEBRTC_SPL_MUL_16_16(*xptr1, *xptr2);
     97                 xptr1++;
     98                 xptr2++;
     99             }
    100 
    101             for (j = loops4; j < loops; j++)
    102             {
    103                 sum += WEBRTC_SPL_MUL_16_16(*xptr1, *xptr2);
    104                 xptr1++;
    105                 xptr2++;
    106             }
    107         }
    108         else
    109         {
    110             for (j = 0; j < loops4; j = j + 4)
    111             {
    112                 sum += WEBRTC_SPL_MUL_16_16_RSFT(*xptr1, *xptr2, scaling);
    113                 xptr1++;
    114                 xptr2++;
    115                 sum += WEBRTC_SPL_MUL_16_16_RSFT(*xptr1, *xptr2, scaling);
    116                 xptr1++;
    117                 xptr2++;
    118                 sum += WEBRTC_SPL_MUL_16_16_RSFT(*xptr1, *xptr2, scaling);
    119                 xptr1++;
    120                 xptr2++;
    121                 sum += WEBRTC_SPL_MUL_16_16_RSFT(*xptr1, *xptr2, scaling);
    122                 xptr1++;
    123                 xptr2++;
    124             }
    125 
    126             for (j = loops4; j < loops; j++)
    127             {
    128                 sum += WEBRTC_SPL_MUL_16_16_RSFT(*xptr1, *xptr2, scaling);
    129                 xptr1++;
    130                 xptr2++;
    131             }
    132         }
    133 
    134 #endif
    135         *resultptr++ = sum;
    136     }
    137 
    138     *scale = scaling;
    139 
    140     return order + 1;
    141 }
    142