Home | History | Annotate | Download | only in libFLAC
      1 /* libFLAC - Free Lossless Audio Codec library
      2  * Copyright (C) 2000-2009  Josh Coalson
      3  * Copyright (C) 2011-2016  Xiph.Org Foundation
      4  *
      5  * Redistribution and use in source and binary forms, with or without
      6  * modification, are permitted provided that the following conditions
      7  * are met:
      8  *
      9  * - Redistributions of source code must retain the above copyright
     10  * notice, this list of conditions and the following disclaimer.
     11  *
     12  * - Redistributions in binary form must reproduce the above copyright
     13  * notice, this list of conditions and the following disclaimer in the
     14  * documentation and/or other materials provided with the distribution.
     15  *
     16  * - Neither the name of the Xiph.org Foundation nor the names of its
     17  * contributors may be used to endorse or promote products derived from
     18  * this software without specific prior written permission.
     19  *
     20  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
     21  * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
     22  * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
     23  * A PARTICULAR PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE FOUNDATION OR
     24  * CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
     25  * EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
     26  * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
     27  * PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
     28  * LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
     29  * NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
     30  * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
     31  */
     32 
     33 #ifdef HAVE_CONFIG_H
     34 #  include <config.h>
     35 #endif
     36 
     37 #include "private/cpu.h"
     38 
     39 #ifndef FLAC__INTEGER_ONLY_LIBRARY
     40 #ifndef FLAC__NO_ASM
     41 #if (defined FLAC__CPU_IA32 || defined FLAC__CPU_X86_64) && defined FLAC__HAS_X86INTRIN
     42 #include "private/fixed.h"
     43 #ifdef FLAC__SSE2_SUPPORTED
     44 
     45 #include <emmintrin.h> /* SSE2 */
     46 #include <math.h>
     47 #include "private/macros.h"
     48 #include "share/compat.h"
     49 #include "FLAC/assert.h"
     50 
     51 #ifdef FLAC__CPU_IA32
     52 #define m128i_to_i64(dest, src) _mm_storel_epi64((__m128i*)&dest, src)
     53 #else
     54 #define m128i_to_i64(dest, src) dest = _mm_cvtsi128_si64(src)
     55 #endif
     56 
     57 FLAC__SSE_TARGET("sse2")
     58 unsigned FLAC__fixed_compute_best_predictor_intrin_sse2(const FLAC__int32 data[], unsigned data_len, float residual_bits_per_sample[FLAC__MAX_FIXED_ORDER + 1])
     59 {
     60 	FLAC__uint32 total_error_0, total_error_1, total_error_2, total_error_3, total_error_4;
     61 	unsigned i, order;
     62 
     63 	__m128i total_err0, total_err1, total_err2;
     64 
     65 	{
     66 		FLAC__int32 itmp;
     67 		__m128i last_error;
     68 
     69 		last_error = _mm_cvtsi32_si128(data[-1]);							// 0   0   0   le0
     70 		itmp = data[-2];
     71 		last_error = _mm_shuffle_epi32(last_error, _MM_SHUFFLE(2,1,0,0));
     72 		last_error = _mm_sub_epi32(last_error, _mm_cvtsi32_si128(itmp));	// 0   0   le0 le1
     73 		itmp -= data[-3];
     74 		last_error = _mm_shuffle_epi32(last_error, _MM_SHUFFLE(2,1,0,0));
     75 		last_error = _mm_sub_epi32(last_error, _mm_cvtsi32_si128(itmp));	// 0   le0 le1 le2
     76 		itmp -= data[-3] - data[-4];
     77 		last_error = _mm_shuffle_epi32(last_error, _MM_SHUFFLE(2,1,0,0));
     78 		last_error = _mm_sub_epi32(last_error, _mm_cvtsi32_si128(itmp));	// le0 le1 le2 le3
     79 
     80 		total_err0 = total_err1 = _mm_setzero_si128();
     81 		for(i = 0; i < data_len; i++) {
     82 			__m128i err0, err1, tmp;
     83 			err0 = _mm_cvtsi32_si128(data[i]);								// 0   0   0   e0
     84 			err1 = _mm_shuffle_epi32(err0, _MM_SHUFFLE(0,0,0,0));			// e0  e0  e0  e0
     85 #if 1 /* OPT_SSE */
     86 			err1 = _mm_sub_epi32(err1, last_error);
     87 			last_error = _mm_srli_si128(last_error, 4);						// 0   le0 le1 le2
     88 			err1 = _mm_sub_epi32(err1, last_error);
     89 			last_error = _mm_srli_si128(last_error, 4);						// 0   0   le0 le1
     90 			err1 = _mm_sub_epi32(err1, last_error);
     91 			last_error = _mm_srli_si128(last_error, 4);						// 0   0   0   le0
     92 			err1 = _mm_sub_epi32(err1, last_error);							// e1  e2  e3  e4
     93 #else
     94 			last_error = _mm_add_epi32(last_error, _mm_srli_si128(last_error, 8));	// le0  le1  le2+le0  le3+le1
     95 			last_error = _mm_add_epi32(last_error, _mm_srli_si128(last_error, 4));	// le0  le1+le0  le2+le0+le1  le3+le1+le2+le0
     96 			err1 = _mm_sub_epi32(err1, last_error);							// e1  e2  e3  e4
     97 #endif
     98 			tmp = _mm_slli_si128(err0, 12);									// e0   0   0   0
     99 			last_error = _mm_srli_si128(err1, 4);							//  0  e1  e2  e3
    100 			last_error = _mm_or_si128(last_error, tmp);						// e0  e1  e2  e3
    101 
    102 			tmp = _mm_srai_epi32(err0, 31);
    103 			err0 = _mm_xor_si128(err0, tmp);
    104 			err0 = _mm_sub_epi32(err0, tmp);
    105 			tmp = _mm_srai_epi32(err1, 31);
    106 			err1 = _mm_xor_si128(err1, tmp);
    107 			err1 = _mm_sub_epi32(err1, tmp);
    108 
    109 			total_err0 = _mm_add_epi32(total_err0, err0);					// 0   0   0   te0
    110 			total_err1 = _mm_add_epi32(total_err1, err1);					// te1 te2 te3 te4
    111 		}
    112 	}
    113 
    114 	total_error_0 = _mm_cvtsi128_si32(total_err0);
    115 	total_err2 = total_err1;											// te1  te2  te3  te4
    116 	total_err1 = _mm_srli_si128(total_err1, 8);							//  0    0   te1  te2
    117 	total_error_4 = _mm_cvtsi128_si32(total_err2);
    118 	total_error_2 = _mm_cvtsi128_si32(total_err1);
    119 	total_err2 = _mm_srli_si128(total_err2,	4);							//  0   te1  te2  te3
    120 	total_err1 = _mm_srli_si128(total_err1, 4);							//  0    0    0   te1
    121 	total_error_3 = _mm_cvtsi128_si32(total_err2);
    122 	total_error_1 = _mm_cvtsi128_si32(total_err1);
    123 
    124 	/* prefer higher order */
    125 	if(total_error_0 < flac_min(flac_min(flac_min(total_error_1, total_error_2), total_error_3), total_error_4))
    126 		order = 0;
    127 	else if(total_error_1 < flac_min(flac_min(total_error_2, total_error_3), total_error_4))
    128 		order = 1;
    129 	else if(total_error_2 < flac_min(total_error_3, total_error_4))
    130 		order = 2;
    131 	else if(total_error_3 < total_error_4)
    132 		order = 3;
    133 	else
    134 		order = 4;
    135 
    136 	/* Estimate the expected number of bits per residual signal sample. */
    137 	/* 'total_error*' is linearly related to the variance of the residual */
    138 	/* signal, so we use it directly to compute E(|x|) */
    139 	FLAC__ASSERT(data_len > 0 || total_error_0 == 0);
    140 	FLAC__ASSERT(data_len > 0 || total_error_1 == 0);
    141 	FLAC__ASSERT(data_len > 0 || total_error_2 == 0);
    142 	FLAC__ASSERT(data_len > 0 || total_error_3 == 0);
    143 	FLAC__ASSERT(data_len > 0 || total_error_4 == 0);
    144 
    145 	residual_bits_per_sample[0] = (float)((total_error_0 > 0) ? log(M_LN2 * (double)total_error_0 / (double)data_len) / M_LN2 : 0.0);
    146 	residual_bits_per_sample[1] = (float)((total_error_1 > 0) ? log(M_LN2 * (double)total_error_1 / (double)data_len) / M_LN2 : 0.0);
    147 	residual_bits_per_sample[2] = (float)((total_error_2 > 0) ? log(M_LN2 * (double)total_error_2 / (double)data_len) / M_LN2 : 0.0);
    148 	residual_bits_per_sample[3] = (float)((total_error_3 > 0) ? log(M_LN2 * (double)total_error_3 / (double)data_len) / M_LN2 : 0.0);
    149 	residual_bits_per_sample[4] = (float)((total_error_4 > 0) ? log(M_LN2 * (double)total_error_4 / (double)data_len) / M_LN2 : 0.0);
    150 
    151 	return order;
    152 }
    153 
    154 FLAC__SSE_TARGET("sse2")
    155 unsigned FLAC__fixed_compute_best_predictor_wide_intrin_sse2(const FLAC__int32 data[], unsigned data_len, float residual_bits_per_sample[FLAC__MAX_FIXED_ORDER + 1])
    156 {
    157 	FLAC__uint64 total_error_0, total_error_1, total_error_2, total_error_3, total_error_4;
    158 	unsigned i, order;
    159 
    160 	__m128i total_err0, total_err1, total_err3;
    161 
    162 	{
    163 		FLAC__int32 itmp;
    164 		__m128i last_error, zero = _mm_setzero_si128();
    165 
    166 		last_error = _mm_cvtsi32_si128(data[-1]);							// 0   0   0   le0
    167 		itmp = data[-2];
    168 		last_error = _mm_shuffle_epi32(last_error, _MM_SHUFFLE(2,1,0,0));
    169 		last_error = _mm_sub_epi32(last_error, _mm_cvtsi32_si128(itmp));	// 0   0   le0 le1
    170 		itmp -= data[-3];
    171 		last_error = _mm_shuffle_epi32(last_error, _MM_SHUFFLE(2,1,0,0));
    172 		last_error = _mm_sub_epi32(last_error, _mm_cvtsi32_si128(itmp));	// 0   le0 le1 le2
    173 		itmp -= data[-3] - data[-4];
    174 		last_error = _mm_shuffle_epi32(last_error, _MM_SHUFFLE(2,1,0,0));
    175 		last_error = _mm_sub_epi32(last_error, _mm_cvtsi32_si128(itmp));	// le0 le1 le2 le3
    176 
    177 		total_err0 = total_err1 = total_err3 = _mm_setzero_si128();
    178 		for(i = 0; i < data_len; i++) {
    179 			__m128i err0, err1, tmp;
    180 			err0 = _mm_cvtsi32_si128(data[i]);								// 0   0   0   e0
    181 			err1 = _mm_shuffle_epi32(err0, _MM_SHUFFLE(0,0,0,0));			// e0  e0  e0  e0
    182 #if 1 /* OPT_SSE */
    183 			err1 = _mm_sub_epi32(err1, last_error);
    184 			last_error = _mm_srli_si128(last_error, 4);						// 0   le0 le1 le2
    185 			err1 = _mm_sub_epi32(err1, last_error);
    186 			last_error = _mm_srli_si128(last_error, 4);						// 0   0   le0 le1
    187 			err1 = _mm_sub_epi32(err1, last_error);
    188 			last_error = _mm_srli_si128(last_error, 4);						// 0   0   0   le0
    189 			err1 = _mm_sub_epi32(err1, last_error);							// e1  e2  e3  e4
    190 #else
    191 			last_error = _mm_add_epi32(last_error, _mm_srli_si128(last_error, 8));	// le0  le1  le2+le0  le3+le1
    192 			last_error = _mm_add_epi32(last_error, _mm_srli_si128(last_error, 4));	// le0  le1+le0  le2+le0+le1  le3+le1+le2+le0
    193 			err1 = _mm_sub_epi32(err1, last_error);							// e1  e2  e3  e4
    194 #endif
    195 			tmp = _mm_slli_si128(err0, 12);									// e0   0   0   0
    196 			last_error = _mm_srli_si128(err1, 4);							//  0  e1  e2  e3
    197 			last_error = _mm_or_si128(last_error, tmp);						// e0  e1  e2  e3
    198 
    199 			tmp = _mm_srai_epi32(err0, 31);
    200 			err0 = _mm_xor_si128(err0, tmp);
    201 			err0 = _mm_sub_epi32(err0, tmp);
    202 			tmp = _mm_srai_epi32(err1, 31);
    203 			err1 = _mm_xor_si128(err1, tmp);
    204 			err1 = _mm_sub_epi32(err1, tmp);
    205 
    206 			total_err0 = _mm_add_epi64(total_err0, err0);					//        0       te0
    207 			err0 = _mm_unpacklo_epi32(err1, zero);							//   0  |e3|   0  |e4|
    208 			err1 = _mm_unpackhi_epi32(err1, zero);							//   0  |e1|   0  |e2|
    209 			total_err3 = _mm_add_epi64(total_err3, err0);					//       te3      te4
    210 			total_err1 = _mm_add_epi64(total_err1, err1);					//       te1      te2
    211 		}
    212 	}
    213 
    214 	m128i_to_i64(total_error_0, total_err0);
    215 	m128i_to_i64(total_error_4, total_err3);
    216 	m128i_to_i64(total_error_2, total_err1);
    217 	total_err3 = _mm_srli_si128(total_err3,	8);							//         0      te3
    218 	total_err1 = _mm_srli_si128(total_err1, 8);							//         0      te1
    219 	m128i_to_i64(total_error_3, total_err3);
    220 	m128i_to_i64(total_error_1, total_err1);
    221 
    222 	/* prefer higher order */
    223 	if(total_error_0 < flac_min(flac_min(flac_min(total_error_1, total_error_2), total_error_3), total_error_4))
    224 		order = 0;
    225 	else if(total_error_1 < flac_min(flac_min(total_error_2, total_error_3), total_error_4))
    226 		order = 1;
    227 	else if(total_error_2 < flac_min(total_error_3, total_error_4))
    228 		order = 2;
    229 	else if(total_error_3 < total_error_4)
    230 		order = 3;
    231 	else
    232 		order = 4;
    233 
    234 	/* Estimate the expected number of bits per residual signal sample. */
    235 	/* 'total_error*' is linearly related to the variance of the residual */
    236 	/* signal, so we use it directly to compute E(|x|) */
    237 	FLAC__ASSERT(data_len > 0 || total_error_0 == 0);
    238 	FLAC__ASSERT(data_len > 0 || total_error_1 == 0);
    239 	FLAC__ASSERT(data_len > 0 || total_error_2 == 0);
    240 	FLAC__ASSERT(data_len > 0 || total_error_3 == 0);
    241 	FLAC__ASSERT(data_len > 0 || total_error_4 == 0);
    242 
    243 	residual_bits_per_sample[0] = (float)((total_error_0 > 0) ? log(M_LN2 * (double)total_error_0 / (double)data_len) / M_LN2 : 0.0);
    244 	residual_bits_per_sample[1] = (float)((total_error_1 > 0) ? log(M_LN2 * (double)total_error_1 / (double)data_len) / M_LN2 : 0.0);
    245 	residual_bits_per_sample[2] = (float)((total_error_2 > 0) ? log(M_LN2 * (double)total_error_2 / (double)data_len) / M_LN2 : 0.0);
    246 	residual_bits_per_sample[3] = (float)((total_error_3 > 0) ? log(M_LN2 * (double)total_error_3 / (double)data_len) / M_LN2 : 0.0);
    247 	residual_bits_per_sample[4] = (float)((total_error_4 > 0) ? log(M_LN2 * (double)total_error_4 / (double)data_len) / M_LN2 : 0.0);
    248 
    249 	return order;
    250 }
    251 
    252 #endif /* FLAC__SSE2_SUPPORTED */
    253 #endif /* (FLAC__CPU_IA32 || FLAC__CPU_X86_64) && FLAC__HAS_X86INTRIN */
    254 #endif /* FLAC__NO_ASM */
    255 #endif /* FLAC__INTEGER_ONLY_LIBRARY */
    256