michael@0: /* Copyright (c) 2013 Jean-Marc Valin and John Ridges */
michael@0: /**
michael@0:    @file pitch_sse.h
michael@0:    @brief Pitch analysis
michael@0:  */
michael@0: 
michael@0: /*
michael@0:    Redistribution and use in source and binary forms, with or without
michael@0:    modification, are permitted provided that the following conditions
michael@0:    are met:
michael@0: 
michael@0:    - Redistributions of source code must retain the above copyright
michael@0:    notice, this list of conditions and the following disclaimer.
michael@0: 
michael@0:    - Redistributions in binary form must reproduce the above copyright
michael@0:    notice, this list of conditions and the following disclaimer in the
michael@0:    documentation and/or other materials provided with the distribution.
michael@0: 
michael@0:    THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
michael@0:    ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
michael@0:    LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
michael@0:    A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER
michael@0:    OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
michael@0:    EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
michael@0:    PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
michael@0:    PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
michael@0:    LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
michael@0:    NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
michael@0:    SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
michael@0: */
michael@0: 
michael@0: #ifndef PITCH_SSE_H
michael@0: #define PITCH_SSE_H
michael@0: 
michael@0: #include <xmmintrin.h>
michael@0: #include "arch.h"
michael@0: 
michael@0: #define OVERRIDE_XCORR_KERNEL
michael@0: static OPUS_INLINE void xcorr_kernel(const opus_val16 *x, const opus_val16 *y, opus_val32 sum[4], int len)
michael@0: {
michael@0:    int j;
michael@0:    __m128 xsum1, xsum2;
michael@0:    xsum1 = _mm_loadu_ps(sum);
michael@0:    xsum2 = _mm_setzero_ps();
michael@0: 
michael@0:    for (j = 0; j < len-3; j += 4)
michael@0:    {
michael@0:       __m128 x0 = _mm_loadu_ps(x+j);
michael@0:       __m128 yj = _mm_loadu_ps(y+j);
michael@0:       __m128 y3 = _mm_loadu_ps(y+j+3);
michael@0: 
michael@0:       xsum1 = _mm_add_ps(xsum1,_mm_mul_ps(_mm_shuffle_ps(x0,x0,0x00),yj));
michael@0:       xsum2 = _mm_add_ps(xsum2,_mm_mul_ps(_mm_shuffle_ps(x0,x0,0x55),
michael@0:                                           _mm_shuffle_ps(yj,y3,0x49)));
michael@0:       xsum1 = _mm_add_ps(xsum1,_mm_mul_ps(_mm_shuffle_ps(x0,x0,0xaa),
michael@0:                                           _mm_shuffle_ps(yj,y3,0x9e)));
michael@0:       xsum2 = _mm_add_ps(xsum2,_mm_mul_ps(_mm_shuffle_ps(x0,x0,0xff),y3));
michael@0:    }
michael@0:    if (j < len)
michael@0:    {
michael@0:       xsum1 = _mm_add_ps(xsum1,_mm_mul_ps(_mm_load1_ps(x+j),_mm_loadu_ps(y+j)));
michael@0:       if (++j < len)
michael@0:       {
michael@0:          xsum2 = _mm_add_ps(xsum2,_mm_mul_ps(_mm_load1_ps(x+j),_mm_loadu_ps(y+j)));
michael@0:          if (++j < len)
michael@0:          {
michael@0:             xsum1 = _mm_add_ps(xsum1,_mm_mul_ps(_mm_load1_ps(x+j),_mm_loadu_ps(y+j)));
michael@0:          }
michael@0:       }
michael@0:    }
michael@0:    _mm_storeu_ps(sum,_mm_add_ps(xsum1,xsum2));
michael@0: }
michael@0: 
michael@0: #define OVERRIDE_DUAL_INNER_PROD
michael@0: static OPUS_INLINE void dual_inner_prod(const opus_val16 *x, const opus_val16 *y01, const opus_val16 *y02,
michael@0:       int N, opus_val32 *xy1, opus_val32 *xy2)
michael@0: {
michael@0:    int i;
michael@0:    __m128 xsum1, xsum2;
michael@0:    xsum1 = _mm_setzero_ps();
michael@0:    xsum2 = _mm_setzero_ps();
michael@0:    for (i=0;i<N-3;i+=4)
michael@0:    {
michael@0:       __m128 xi = _mm_loadu_ps(x+i);
michael@0:       __m128 y1i = _mm_loadu_ps(y01+i);
michael@0:       __m128 y2i = _mm_loadu_ps(y02+i);
michael@0:       xsum1 = _mm_add_ps(xsum1,_mm_mul_ps(xi, y1i));
michael@0:       xsum2 = _mm_add_ps(xsum2,_mm_mul_ps(xi, y2i));
michael@0:    }
michael@0:    /* Horizontal sum */
michael@0:    xsum1 = _mm_add_ps(xsum1, _mm_movehl_ps(xsum1, xsum1));
michael@0:    xsum1 = _mm_add_ss(xsum1, _mm_shuffle_ps(xsum1, xsum1, 0x55));
michael@0:    _mm_store_ss(xy1, xsum1);
michael@0:    xsum2 = _mm_add_ps(xsum2, _mm_movehl_ps(xsum2, xsum2));
michael@0:    xsum2 = _mm_add_ss(xsum2, _mm_shuffle_ps(xsum2, xsum2, 0x55));
michael@0:    _mm_store_ss(xy2, xsum2);
michael@0:    for (;i<N;i++)
michael@0:    {
michael@0:       *xy1 = MAC16_16(*xy1, x[i], y01[i]);
michael@0:       *xy2 = MAC16_16(*xy2, x[i], y02[i]);
michael@0:    }
michael@0: }
michael@0: 
michael@0: #define OVERRIDE_COMB_FILTER_CONST
michael@0: static OPUS_INLINE void comb_filter_const(opus_val32 *y, opus_val32 *x, int T, int N,
michael@0:       opus_val16 g10, opus_val16 g11, opus_val16 g12)
michael@0: {
michael@0:    int i;
michael@0:    __m128 x0v;
michael@0:    __m128 g10v, g11v, g12v;
michael@0:    g10v = _mm_load1_ps(&g10);
michael@0:    g11v = _mm_load1_ps(&g11);
michael@0:    g12v = _mm_load1_ps(&g12);
michael@0:    x0v = _mm_loadu_ps(&x[-T-2]);
michael@0:    for (i=0;i<N-3;i+=4)
michael@0:    {
michael@0:       __m128 yi, yi2, x1v, x2v, x3v, x4v;
michael@0:       const opus_val32 *xp = &x[i-T-2];
michael@0:       yi = _mm_loadu_ps(x+i);
michael@0:       x4v = _mm_loadu_ps(xp+4);
michael@0: #if 0
michael@0:       /* Slower version with all loads */
michael@0:       x1v = _mm_loadu_ps(xp+1);
michael@0:       x2v = _mm_loadu_ps(xp+2);
michael@0:       x3v = _mm_loadu_ps(xp+3);
michael@0: #else
michael@0:       x2v = _mm_shuffle_ps(x0v, x4v, 0x4e);
michael@0:       x1v = _mm_shuffle_ps(x0v, x2v, 0x99);
michael@0:       x3v = _mm_shuffle_ps(x2v, x4v, 0x99);
michael@0: #endif
michael@0: 
michael@0:       yi = _mm_add_ps(yi, _mm_mul_ps(g10v,x2v));
michael@0: #if 0 /* Set to 1 to make it bit-exact with the non-SSE version */
michael@0:       yi = _mm_add_ps(yi, _mm_mul_ps(g11v,_mm_add_ps(x3v,x1v)));
michael@0:       yi = _mm_add_ps(yi, _mm_mul_ps(g12v,_mm_add_ps(x4v,x0v)));
michael@0: #else
michael@0:       /* Use partial sums */
michael@0:       yi2 = _mm_add_ps(_mm_mul_ps(g11v,_mm_add_ps(x3v,x1v)),
michael@0:                        _mm_mul_ps(g12v,_mm_add_ps(x4v,x0v)));
michael@0:       yi = _mm_add_ps(yi, yi2);
michael@0: #endif
michael@0:       x0v=x4v;
michael@0:       _mm_storeu_ps(y+i, yi);
michael@0:    }
michael@0: #ifdef CUSTOM_MODES
michael@0:    for (;i<N;i++)
michael@0:    {
michael@0:       y[i] = x[i]
michael@0:                + MULT16_32_Q15(g10,x[i-T])
michael@0:                + MULT16_32_Q15(g11,ADD32(x[i-T+1],x[i-T-1]))
michael@0:                + MULT16_32_Q15(g12,ADD32(x[i-T+2],x[i-T-2]));
michael@0:    }
michael@0: #endif
michael@0: }
michael@0: 
michael@0: #endif