Upgrade openjpeg to revision 2997.

core/src/fxcodec/fx_libopenjpeg/libopenjpeg20/mct.c

 /*
  * The copyright in this software is being made available under the 2-clauses 
  * BSD License, included below. This software may be subject to other third 
  * party and contributor rights, including patent rights, and no such rights
  * are granted under this license.
  *
  * Copyright (c) 2002-2014, Universite catholique de Louvain (UCL), Belgium
  * Copyright (c) 2002-2014, Professor Benoit Macq
  * Copyright (c) 2001-2003, David Janssens
  * Copyright (c) 2002-2003, Yannick Verschueren
@@ skipping 22 matching lines @@
  * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
  * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
  * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
  * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
  * POSSIBILITY OF SUCH DAMAGE.
  */
 
 #ifdef __SSE__
 #include <xmmintrin.h>
 #endif
+#ifdef __SSE2__
+#include <emmintrin.h>
+#endif
+#ifdef __SSE4_1__
+#include <smmintrin.h>
+#endif
 
 #include "opj_includes.h"
 
 /* <summary> */
 /* This table contains the norms of the basis function of the reversible MCT. */
 /* </summary> */
 static const OPJ_FLOAT64 opj_mct_norms[3] = { 1.732, .8292, .8292 };
 
 /* <summary> */
 /* This table contains the norms of the basis function of the irreversible MCT. */
 /* </summary> */
 static const OPJ_FLOAT64 opj_mct_norms_real[3] = { 1.732, 1.805, 1.573 };
 
 const OPJ_FLOAT64 * opj_mct_get_mct_norms ()
 {
         return opj_mct_norms;
 }
 
 const OPJ_FLOAT64 * opj_mct_get_mct_norms_real ()
 {
         return opj_mct_norms_real;
 }
 
 /* <summary> */
 /* Foward reversible MCT. */
 /* </summary> */
+#ifdef __SSE2__
 void opj_mct_encode(
                 OPJ_INT32* restrict c0,
                 OPJ_INT32* restrict c1,
                 OPJ_INT32* restrict c2,
                 OPJ_UINT32 n)
 {
-        OPJ_UINT32 i;
-        for(i = 0; i < n; ++i) {
+        OPJ_SIZE_T i;
+        const OPJ_SIZE_T len = n;
+        
+        for(i = 0; i < (len & ~3U); i += 4) {
+                __m128i y, u, v;
+                __m128i r = _mm_load_si128((const __m128i *)&(c0[i]));
+                __m128i g = _mm_load_si128((const __m128i *)&(c1[i]));
+                __m128i b = _mm_load_si128((const __m128i *)&(c2[i]));
+                y = _mm_add_epi32(g, g);
+                y = _mm_add_epi32(y, b);
+                y = _mm_add_epi32(y, r);
+                y = _mm_srai_epi32(y, 2);
+                u = _mm_sub_epi32(b, g);
+                v = _mm_sub_epi32(r, g);
+                _mm_store_si128((__m128i *)&(c0[i]), y);
+                _mm_store_si128((__m128i *)&(c1[i]), u);
+                _mm_store_si128((__m128i *)&(c2[i]), v);
+        }
+        
+        for(; i < len; ++i) {
                 OPJ_INT32 r = c0[i];
                 OPJ_INT32 g = c1[i];
                 OPJ_INT32 b = c2[i];
                 OPJ_INT32 y = (r + (g * 2) + b) >> 2;
                 OPJ_INT32 u = b - g;
                 OPJ_INT32 v = r - g;
                 c0[i] = y;
                 c1[i] = u;
                 c2[i] = v;
         }
 }
+#else
+void opj_mct_encode(
+                OPJ_INT32* restrict c0,
+                OPJ_INT32* restrict c1,
+                OPJ_INT32* restrict c2,
+                OPJ_UINT32 n)
+{
+        OPJ_SIZE_T i;
+        const OPJ_SIZE_T len = n;
+        
+        for(i = 0; i < len; ++i) {
+                OPJ_INT32 r = c0[i];
+                OPJ_INT32 g = c1[i];
+                OPJ_INT32 b = c2[i];
+                OPJ_INT32 y = (r + (g * 2) + b) >> 2;
+                OPJ_INT32 u = b - g;
+                OPJ_INT32 v = r - g;
+                c0[i] = y;
+                c1[i] = u;
+                c2[i] = v;
+        }
+}
+#endif
 
 /* <summary> */
 /* Inverse reversible MCT. */
 /* </summary> */
+#ifdef __SSE2__
+void opj_mct_decode(
+                OPJ_INT32* restrict c0,
+                OPJ_INT32* restrict c1,
+                OPJ_INT32* restrict c2,
+                OPJ_UINT32 n)
+{
+        OPJ_SIZE_T i;
+        const OPJ_SIZE_T len = n;
+        
+        for(i = 0; i < (len & ~3U); i += 4) {
+                __m128i r, g, b;
+                __m128i y = _mm_load_si128((const __m128i *)&(c0[i]));
+                __m128i u = _mm_load_si128((const __m128i *)&(c1[i]));
+                __m128i v = _mm_load_si128((const __m128i *)&(c2[i]));
+                g = y;
+                g = _mm_sub_epi32(g, _mm_srai_epi32(_mm_add_epi32(u, v), 2));
+                r = _mm_add_epi32(v, g);
+                b = _mm_add_epi32(u, g);
+                _mm_store_si128((__m128i *)&(c0[i]), r);
+                _mm_store_si128((__m128i *)&(c1[i]), g);
+                _mm_store_si128((__m128i *)&(c2[i]), b);
+        }
+        for (; i < len; ++i) {
+                OPJ_INT32 y = c0[i];
+                OPJ_INT32 u = c1[i];
+                OPJ_INT32 v = c2[i];
+                OPJ_INT32 g = y - ((u + v) >> 2);
+                OPJ_INT32 r = v + g;
+                OPJ_INT32 b = u + g;
+                c0[i] = r;
+                c1[i] = g;
+                c2[i] = b;
+        }
+}
+#else
 void opj_mct_decode(
                 OPJ_INT32* restrict c0,
                 OPJ_INT32* restrict c1, 
                 OPJ_INT32* restrict c2, 
                 OPJ_UINT32 n)
 {
         OPJ_UINT32 i;
         for (i = 0; i < n; ++i) {
                 OPJ_INT32 y = c0[i];
                 OPJ_INT32 u = c1[i];
                 OPJ_INT32 v = c2[i];
                 OPJ_INT32 g = y - ((u + v) >> 2);
                 OPJ_INT32 r = v + g;
                 OPJ_INT32 b = u + g;
                 c0[i] = r;
                 c1[i] = g;
                 c2[i] = b;
         }
 }
+#endif
 
 /* <summary> */
 /* Get norm of basis function of reversible MCT. */
 /* </summary> */
 OPJ_FLOAT64 opj_mct_getnorm(OPJ_UINT32 compno) {
         return opj_mct_norms[compno];
 }
 
 /* <summary> */
 /* Foward irreversible MCT. */
 /* </summary> */
+#ifdef __SSE4_1__
+void opj_mct_encode_real(
+                                                                                                 OPJ_INT32* restrict c0,
+                                                                                                 OPJ_INT32* restrict c1,
+                                                                                                 OPJ_INT32* restrict c2,
+                                                                                                 OPJ_UINT32 n)
+{
+        OPJ_SIZE_T i;
+        const OPJ_SIZE_T len = n;
+        
+        const __m128i ry = _mm_set1_epi32(2449);
+        const __m128i gy = _mm_set1_epi32(4809);
+        const __m128i by = _mm_set1_epi32(934);
+        const __m128i ru = _mm_set1_epi32(1382);
+        const __m128i gu = _mm_set1_epi32(2714);
+        /* const __m128i bu = _mm_set1_epi32(4096); */
+        /* const __m128i rv = _mm_set1_epi32(4096); */
+        const __m128i gv = _mm_set1_epi32(3430);
+        const __m128i bv = _mm_set1_epi32(666);
+        const __m128i mulround = _mm_shuffle_epi32(_mm_cvtsi32_si128(4096), _MM_SHUFFLE(1, 0, 1, 0));
+        
+        for(i = 0; i < (len & ~3U); i += 4) {
+                __m128i lo, hi;
+                __m128i y, u, v;
+                __m128i r = _mm_load_si128((const __m128i *)&(c0[i]));
+                __m128i g = _mm_load_si128((const __m128i *)&(c1[i]));
+                __m128i b = _mm_load_si128((const __m128i *)&(c2[i]));
+                
+                lo = r;
+                hi = _mm_shuffle_epi32(r, _MM_SHUFFLE(3, 3, 1, 1));
+                lo = _mm_mul_epi32(lo, ry);
+                hi = _mm_mul_epi32(hi, ry);
+                lo = _mm_add_epi64(lo, mulround);
+                hi = _mm_add_epi64(hi, mulround);
+                lo = _mm_srli_epi64(lo, 13);
+                hi = _mm_slli_epi64(hi, 32-13);
+                y = _mm_blend_epi16(lo, hi, 0xCC);
+                
+                lo = g;
+                hi = _mm_shuffle_epi32(g, _MM_SHUFFLE(3, 3, 1, 1));
+                lo = _mm_mul_epi32(lo, gy);
+                hi = _mm_mul_epi32(hi, gy);
+                lo = _mm_add_epi64(lo, mulround);
+                hi = _mm_add_epi64(hi, mulround);
+                lo = _mm_srli_epi64(lo, 13);
+                hi = _mm_slli_epi64(hi, 32-13);
+                y = _mm_add_epi32(y, _mm_blend_epi16(lo, hi, 0xCC));
+                
+                lo = b;
+                hi = _mm_shuffle_epi32(b, _MM_SHUFFLE(3, 3, 1, 1));
+                lo = _mm_mul_epi32(lo, by);
+                hi = _mm_mul_epi32(hi, by);
+                lo = _mm_add_epi64(lo, mulround);
+                hi = _mm_add_epi64(hi, mulround);
+                lo = _mm_srli_epi64(lo, 13);
+                hi = _mm_slli_epi64(hi, 32-13);
+                y = _mm_add_epi32(y, _mm_blend_epi16(lo, hi, 0xCC));
+                _mm_store_si128((__m128i *)&(c0[i]), y);
+                
+                /*lo = b;
+                hi = _mm_shuffle_epi32(b, _MM_SHUFFLE(3, 3, 1, 1));
+                lo = _mm_mul_epi32(lo, mulround);
+                hi = _mm_mul_epi32(hi, mulround);*/
+                lo = _mm_cvtepi32_epi64(_mm_shuffle_epi32(b, _MM_SHUFFLE(3, 2, 2, 0)));
+                hi = _mm_cvtepi32_epi64(_mm_shuffle_epi32(b, _MM_SHUFFLE(3, 2, 3, 1)));
+                lo = _mm_slli_epi64(lo, 12);
+                hi = _mm_slli_epi64(hi, 12);
+                lo = _mm_add_epi64(lo, mulround);
+                hi = _mm_add_epi64(hi, mulround);
+                lo = _mm_srli_epi64(lo, 13);
+                hi = _mm_slli_epi64(hi, 32-13);
+                u = _mm_blend_epi16(lo, hi, 0xCC);
+                
+                lo = r;
+                hi = _mm_shuffle_epi32(r, _MM_SHUFFLE(3, 3, 1, 1));
+                lo = _mm_mul_epi32(lo, ru);
+                hi = _mm_mul_epi32(hi, ru);
+                lo = _mm_add_epi64(lo, mulround);
+                hi = _mm_add_epi64(hi, mulround);
+                lo = _mm_srli_epi64(lo, 13);
+                hi = _mm_slli_epi64(hi, 32-13);
+                u = _mm_sub_epi32(u, _mm_blend_epi16(lo, hi, 0xCC));
+                
+                lo = g;
+                hi = _mm_shuffle_epi32(g, _MM_SHUFFLE(3, 3, 1, 1));
+                lo = _mm_mul_epi32(lo, gu);
+                hi = _mm_mul_epi32(hi, gu);
+                lo = _mm_add_epi64(lo, mulround);
+                hi = _mm_add_epi64(hi, mulround);
+                lo = _mm_srli_epi64(lo, 13);
+                hi = _mm_slli_epi64(hi, 32-13);
+                u = _mm_sub_epi32(u, _mm_blend_epi16(lo, hi, 0xCC));
+                _mm_store_si128((__m128i *)&(c1[i]), u);
+                
+                /*lo = r;
+                hi = _mm_shuffle_epi32(r, _MM_SHUFFLE(3, 3, 1, 1));
+                lo = _mm_mul_epi32(lo, mulround);
+                hi = _mm_mul_epi32(hi, mulround);*/
+                lo = _mm_cvtepi32_epi64(_mm_shuffle_epi32(r, _MM_SHUFFLE(3, 2, 2, 0)));
+                hi = _mm_cvtepi32_epi64(_mm_shuffle_epi32(r, _MM_SHUFFLE(3, 2, 3, 1)));
+                lo = _mm_slli_epi64(lo, 12);
+                hi = _mm_slli_epi64(hi, 12);
+                lo = _mm_add_epi64(lo, mulround);
+                hi = _mm_add_epi64(hi, mulround);
+                lo = _mm_srli_epi64(lo, 13);
+                hi = _mm_slli_epi64(hi, 32-13);
+                v = _mm_blend_epi16(lo, hi, 0xCC);
+                
+                lo = g;
+                hi = _mm_shuffle_epi32(g, _MM_SHUFFLE(3, 3, 1, 1));
+                lo = _mm_mul_epi32(lo, gv);
+                hi = _mm_mul_epi32(hi, gv);
+                lo = _mm_add_epi64(lo, mulround);
+                hi = _mm_add_epi64(hi, mulround);
+                lo = _mm_srli_epi64(lo, 13);
+                hi = _mm_slli_epi64(hi, 32-13);
+                v = _mm_sub_epi32(v, _mm_blend_epi16(lo, hi, 0xCC));
+                
+                lo = b;
+                hi = _mm_shuffle_epi32(b, _MM_SHUFFLE(3, 3, 1, 1));
+                lo = _mm_mul_epi32(lo, bv);
+                hi = _mm_mul_epi32(hi, bv);
+                lo = _mm_add_epi64(lo, mulround);
+                hi = _mm_add_epi64(hi, mulround);
+                lo = _mm_srli_epi64(lo, 13);
+                hi = _mm_slli_epi64(hi, 32-13);
+                v = _mm_sub_epi32(v, _mm_blend_epi16(lo, hi, 0xCC));
+                _mm_store_si128((__m128i *)&(c2[i]), v);
+        }
+        for(; i < len; ++i) {
+                OPJ_INT32 r = c0[i];
+                OPJ_INT32 g = c1[i];
+                OPJ_INT32 b = c2[i];
+                OPJ_INT32 y =  opj_int_fix_mul(r, 2449) + opj_int_fix_mul(g, 4809) + opj_int_fix_mul(b, 934);
+                OPJ_INT32 u = -opj_int_fix_mul(r, 1382) - opj_int_fix_mul(g, 2714) + opj_int_fix_mul(b, 4096);
+                OPJ_INT32 v =  opj_int_fix_mul(r, 4096) - opj_int_fix_mul(g, 3430) - opj_int_fix_mul(b, 666);
+                c0[i] = y;
+                c1[i] = u;
+                c2[i] = v;
+        }
+}
+#else
 void opj_mct_encode_real(
                 OPJ_INT32* restrict c0,
                 OPJ_INT32* restrict c1,
                 OPJ_INT32* restrict c2,
                 OPJ_UINT32 n)
 {
         OPJ_UINT32 i;
         for(i = 0; i < n; ++i) {
                 OPJ_INT32 r = c0[i];
                 OPJ_INT32 g = c1[i];
                 OPJ_INT32 b = c2[i];
                 OPJ_INT32 y =  opj_int_fix_mul(r, 2449) + opj_int_fix_mul(g, 4809) + opj_int_fix_mul(b, 934);
                 OPJ_INT32 u = -opj_int_fix_mul(r, 1382) - opj_int_fix_mul(g, 2714) + opj_int_fix_mul(b, 4096);
                 OPJ_INT32 v =  opj_int_fix_mul(r, 4096) - opj_int_fix_mul(g, 3430) - opj_int_fix_mul(b, 666);
                 c0[i] = y;
                 c1[i] = u;
                 c2[i] = v;
         }
 }
+#endif
 
 /* <summary> */
 /* Inverse irreversible MCT. */
 /* </summary> */
 void opj_mct_decode_real(
                 OPJ_FLOAT32* restrict c0,
                 OPJ_FLOAT32* restrict c1,
                 OPJ_FLOAT32* restrict c2,
                 OPJ_UINT32 n)
 {
@@ skipping 165 matching lines @@
                 lIndex = i;
 
                 for     (j=0;j<pNbComps;++j) {
                         lCurrentValue = lMatrix[lIndex];
                         lIndex += pNbComps;
                         lNorms[i] += lCurrentValue * lCurrentValue;
                 }
                 lNorms[i] = sqrt(lNorms[i]);
         }
 }