diff --git a/include/libyuv/compare_row.h b/include/libyuv/compare_row.h index 26a62e0ce..22336f9a1 100644 --- a/include/libyuv/compare_row.h +++ b/include/libyuv/compare_row.h @@ -21,49 +21,55 @@ namespace libyuv { extern "C" { #endif -// The following are available for Visual C 32/64 bit, and for clang with -// LIBYUV_ENABLE_ROWWIN, which is also enabled for MemorySanitizer. These are -// intrinsics in compare_win.cc. -#if !defined(LIBYUV_DISABLE_X86) && \ - (defined(__x86_64__) || defined(__i386__) || defined(_M_X64) || \ - defined(_M_IX86)) && \ - ((defined(_MSC_VER) && !defined(__clang__)) || \ - defined(LIBYUV_ENABLE_ROWWIN)) +// The following are available for 32 bit x86, Visual C 64 bit, and for clang +// with LIBYUV_ENABLE_ROWWIN, which is also enabled for MemorySanitizer. These +// are intrinsics in compare_win.cc. +#if !defined(LIBYUV_DISABLE_X86) && \ + (defined(__i386__) || defined(_M_IX86) || \ + ((defined(__x86_64__) || defined(_M_X64)) && \ + ((defined(_MSC_VER) && !defined(__clang__)) || \ + defined(LIBYUV_ENABLE_ROWWIN)))) #define HAS_HAMMINGDISTANCE_SSE42 #define HAS_SUMSQUAREERROR_SSE2 #define HAS_HASHDJB2_SSE41 #define HAS_SSIM4X4ROW_SSE2 +#define HAS_SSIMWINDOWROW_SSE2 #define HAS_HAMMINGDISTANCE_AVX2 #define HAS_SUMSQUAREERROR_AVX2 #define HAS_HASHDJB2_AVX2 #define HAS_SSIM4X4ROW_AVX2 +#define HAS_SSIMWINDOWROW_AVX2 #define HAS_HAMMINGDISTANCE_AVX512BW #define HAS_SUMSQUAREERROR_AVX512BW #define HAS_HASHDJB2_AVX512BW #define HAS_SSIM4X4ROW_AVX512BW +#define HAS_SSIMWINDOWROW_AVX512BW #endif -// The following are available for GCC and clang 32/64 bit, as inline assembly -// in compare_gcc.cc: -#if !defined(LIBYUV_DISABLE_X86) && \ - (defined(__x86_64__) || defined(__i386__)) && \ - !defined(LIBYUV_ENABLE_ROWWIN) +// The following are available for 64 bit GCC and clang, as inline assembly in +// compare_gcc.cc: +#if !defined(LIBYUV_DISABLE_X86) && defined(__x86_64__) && \ + !defined(LIBYUV_ENABLE_ROWWIN) && \ + (!defined(_MSC_VER) || defined(__clang__)) #define HAS_HAMMINGDISTANCE_SSE42 #define HAS_HAMMINGDISTANCE_SSSE3 #define HAS_SUMSQUAREERROR_SSE2 #define HAS_HASHDJB2_SSE41 #define HAS_SSIM4X4ROW_SSE2 +#define HAS_SSIMWINDOWROW_SSE2 #if defined(CLANG_HAS_AVX2) || defined(GCC_HAS_AVX2) #define HAS_HAMMINGDISTANCE_AVX2 #define HAS_SUMSQUAREERROR_AVX2 #define HAS_HASHDJB2_AVX2 #define HAS_SSIM4X4ROW_AVX2 +#define HAS_SSIMWINDOWROW_AVX2 #endif -#if defined(CLANG_HAS_AVX512) +#if defined(CLANG_HAS_AVX512) || defined(GCC_HAS_AVX512) #define HAS_HAMMINGDISTANCE_AVX512BW #define HAS_SUMSQUAREERROR_AVX512BW #define HAS_HASHDJB2_AVX512BW #define HAS_SSIM4X4ROW_AVX512BW +#define HAS_SSIMWINDOWROW_AVX512BW #endif #endif @@ -265,9 +271,18 @@ void Ssim4x4Row_HVX(const uint8_t* src_a, int stride_b, uint32_t* dst, int width); -double SsimWindowRow_C(const uint32_t* row0, - const uint32_t* row1, - int num_windows); +float SsimWindowRow_C(const uint32_t* row0, + const uint32_t* row1, + int num_windows); +float SsimWindowRow_SSE2(const uint32_t* row0, + const uint32_t* row1, + int num_windows); +float SsimWindowRow_AVX2(const uint32_t* row0, + const uint32_t* row1, + int num_windows); +float SsimWindowRow_AVX512BW(const uint32_t* row0, + const uint32_t* row1, + int num_windows); #ifdef __cplusplus } // extern "C" diff --git a/include/libyuv/cpu_support.h b/include/libyuv/cpu_support.h index adb50c42d..1dd9d71a4 100644 --- a/include/libyuv/cpu_support.h +++ b/include/libyuv/cpu_support.h @@ -70,6 +70,15 @@ extern "C" { #endif // clang >= 3.4 #endif // __clang__ +// GCC >= 6.0.0 required for AVX512. +#if defined(__GNUC__) && !defined(__clang__) && \ + !defined(LIBYUV_ENABLE_ROWWIN) && \ + (defined(__x86_64__) || defined(__i386__)) +#if (__GNUC__ >= 6) +#define GCC_HAS_AVX512 1 +#endif // GNUC >= 6 +#endif // __GNUC__ + // clang >= 6.0.0 required for AVX512. #if defined(__clang__) && !defined(LIBYUV_ENABLE_ROWWIN) && \ (defined(__x86_64__) || defined(__i386__)) diff --git a/linux.mk b/linux.mk index 8ab3be195..f4d9c73b1 100644 --- a/linux.mk +++ b/linux.mk @@ -71,11 +71,15 @@ LOCAL_OBJ_FILES := \ .c.o: $(CC) -c $(CFLAGS) $*.c -o $*.o -all: libyuv.a i444tonv12_eg yuvconvert yuvconstants cpuid psnr +all: libyuv.a i444tonv12_eg yuvconvert yuvconstants cpuid psnr compare libyuv.a: $(LOCAL_OBJ_FILES) $(AR) $(ARFLAGS) $@ $(LOCAL_OBJ_FILES) +# A C++ test utility that compares two YUV files. +compare: util/compare.cc libyuv.a + $(CXX) $(CXXFLAGS) -Iutil/ -o $@ util/compare.cc libyuv.a + # A C++ test utility that uses libyuv conversion. yuvconvert: util/yuvconvert.cc libyuv.a $(CXX) $(CXXFLAGS) -Iutil/ -o $@ util/yuvconvert.cc libyuv.a @@ -99,4 +103,4 @@ cpuid: util/cpuid.c libyuv.a $(CC) $(CFLAGS) -o $@ util/cpuid.c libyuv.a clean: - /bin/rm -f source/*.o *.ii *.s libyuv.a i444tonv12_eg yuvconvert yuvconstants cpuid psnr + /bin/rm -f source/*.o *.ii *.s libyuv.a i444tonv12_eg yuvconvert yuvconstants cpuid psnr compare diff --git a/source/compare.cc b/source/compare.cc index 2ec13abc5..0783cdf52 100644 --- a/source/compare.cc +++ b/source/compare.cc @@ -380,18 +380,14 @@ double I420Psnr(const uint8_t* src_y_a, return SumSquareErrorToPsnr(sse, samples); } -static const int64_t cc1 = 26634; // (64^2*(.01*255)^2 -static const int64_t cc2 = 239708; // (64^2*(.03*255)^2 - -static double Ssim8x8_C(const uint8_t* src_a, - int stride_a, - const uint8_t* src_b, - int stride_b) { - int64_t sum_a = 0; - int64_t sum_b = 0; - int64_t sum_sq_a = 0; - int64_t sum_sq_b = 0; - int64_t sum_axb = 0; +static float Ssim8x8_C(const uint8_t* src_a, + int stride_a, + const uint8_t* src_b, + int stride_b) { + int32_t sum_a = 0; + int32_t sum_b = 0; + int32_t sum_sq_ab = 0; + int32_t sum_axb = 0; int i; for (i = 0; i < 8; ++i) { @@ -399,8 +395,7 @@ static double Ssim8x8_C(const uint8_t* src_a, for (j = 0; j < 8; ++j) { sum_a += src_a[j]; sum_b += src_b[j]; - sum_sq_a += src_a[j] * src_a[j]; - sum_sq_b += src_b[j] * src_b[j]; + sum_sq_ab += src_a[j] * src_a[j] + src_b[j] * src_b[j]; sum_axb += src_a[j] * src_b[j]; } @@ -409,27 +404,15 @@ static double Ssim8x8_C(const uint8_t* src_a, } { - const int64_t count = 64; - // scale the constants by number of pixels - const int64_t c1 = (cc1 * count * count) >> 12; - const int64_t c2 = (cc2 * count * count) >> 12; - - const int64_t sum_a_x_sum_b = sum_a * sum_b; - - const int64_t ssim_n = (2 * sum_a_x_sum_b + c1) * - (2 * count * sum_axb - 2 * sum_a_x_sum_b + c2); - - const int64_t sum_a_sq = sum_a * sum_a; - const int64_t sum_b_sq = sum_b * sum_b; - - const int64_t ssim_d = - (sum_a_sq + sum_b_sq + c1) * - (count * sum_sq_a - sum_a_sq + count * sum_sq_b - sum_b_sq + c2); - - if (ssim_d == 0) { - return DBL_MAX; - } - return (double)ssim_n / (double)ssim_d; + const int32_t c1 = 26634; // (64^2*(.01*255)^2) + const int32_t c2 = 239708; // (64^2*(.03*255)^2) + const int32_t sum_a_x_sum_b_2 = 2 * sum_a * sum_b; + const int32_t sum_a_sq_b_sq = sum_a * sum_a + sum_b * sum_b; + const int32_t n1 = sum_a_x_sum_b_2 + c1; + const int32_t n2 = (sum_axb << 7) - sum_a_x_sum_b_2 + c2; + const int32_t d1 = sum_a_sq_b_sq + c1; + const int32_t d2 = (sum_sq_ab << 6) - sum_a_sq_b_sq + c2; + return ((float)n1 * (float)n2) / ((float)d1 * (float)d2); } } @@ -454,6 +437,9 @@ double CalcFrameSsim(const uint8_t* src_a, int simd_size = 4; void (*Ssim4x4Row)(const uint8_t* src_a, int stride_a, const uint8_t* src_b, int stride_b, uint32_t* dst, int width) = Ssim4x4Row_C; + int win_simd_size = 1; + float (*SsimWindowRow)(const uint32_t* row0, const uint32_t* row1, + int num_windows) = SsimWindowRow_C; #if defined(HAS_SSIM4X4ROW_SSE2) if (TestCpuFlag(kCpuHasSSE2)) { Ssim4x4Row = Ssim4x4Row_SSE2; @@ -472,6 +458,24 @@ double CalcFrameSsim(const uint8_t* src_a, simd_size = 4; } #endif +#if defined(HAS_SSIMWINDOWROW_SSE2) + if (TestCpuFlag(kCpuHasSSE2)) { + SsimWindowRow = SsimWindowRow_SSE2; + win_simd_size = 4; + } +#endif +#if defined(HAS_SSIMWINDOWROW_AVX2) + if (TestCpuFlag(kCpuHasAVX2)) { + SsimWindowRow = SsimWindowRow_AVX2; + win_simd_size = 8; + } +#endif +#if defined(HAS_SSIMWINDOWROW_AVX512BW) + if (TestCpuFlag(kCpuHasAVX512BW)) { + SsimWindowRow = SsimWindowRow_AVX512BW; + win_simd_size = 16; + } +#endif #if defined(HAS_SSIM4X4ROW_NEON) if (TestCpuFlag(kCpuHasNEON)) { Ssim4x4Row = Ssim4x4Row_NEON; @@ -532,6 +536,8 @@ double CalcFrameSsim(const uint8_t* src_a, const int tile_width = num_tiles * 4; const int main_width = tile_width & ~(simd_size - 1); const int rem_width = tile_width & (simd_size - 1); + const int main_windows = num_windows & ~(win_simd_size - 1); + const int rem_windows = num_windows & (win_simd_size - 1); #ifdef _OPENMP #pragma omp parallel reduction(+ : ssim_total, samples) @@ -567,7 +573,13 @@ double CalcFrameSsim(const uint8_t* src_a, Ssim4x4Row_C(a + main_width, stride_a, b + main_width, stride_b, row1 + main_width, rem_width); } - ssim_total += SsimWindowRow_C(row0, row1, num_windows); + if (main_windows > 0) { + ssim_total += SsimWindowRow(row0, row1, main_windows); + } + if (rem_windows > 0) { + ssim_total += SsimWindowRow_C(row0 + main_windows * 4, + row1 + main_windows * 4, rem_windows); + } samples += num_windows; uint32_t* tmp = row0; row0 = row1; diff --git a/source/compare_common.cc b/source/compare_common.cc index 5a0533d91..cd667fdd2 100644 --- a/source/compare_common.cc +++ b/source/compare_common.cc @@ -130,46 +130,41 @@ void Ssim4x4Row_C(const uint8_t* src_a, // Combine two consecutive 4x4 tile rows (row0 and row1) into 8x8 windows // stepped on a 4-pixel horizontal grid and return the sum of SSIM values. -double SsimWindowRow_C(const uint32_t* row0, - const uint32_t* row1, - int num_windows) { - const int64_t count = 64; - // scale the constants by number of pixels - const int64_t cc1 = 26634; // (64^2*(.01*255)^2) - const int64_t cc2 = 239708; // (64^2*(.03*255)^2) - const int64_t c1 = (cc1 * count * count) >> 12; - const int64_t c2 = (cc2 * count * count) >> 12; +float SsimWindowRow_C(const uint32_t* row0, + const uint32_t* row1, + int num_windows) { + const int32_t c1 = 26634; // (64^2*(.01*255)^2) + const int32_t c2 = 239708; // (64^2*(.03*255)^2) - double ssim_total = 0.0; - uint32_t prev_a = row0[0] + row1[0]; - uint32_t prev_b = row0[1] + row1[1]; - uint32_t prev_sq = row0[2] + row1[2]; - uint32_t prev_axb = row0[3] + row1[3]; + float ssim_total = 0.0f; + int32_t prev_a = (int32_t)(row0[0] + row1[0]); + int32_t prev_b = (int32_t)(row0[1] + row1[1]); + int32_t prev_sq = (int32_t)(row0[2] + row1[2]); + int32_t prev_axb = (int32_t)(row0[3] + row1[3]); int w; for (w = 0; w < num_windows; ++w) { row0 += 4; row1 += 4; - const uint32_t next_a = row0[0] + row1[0]; - const uint32_t next_b = row0[1] + row1[1]; - const uint32_t next_sq = row0[2] + row1[2]; - const uint32_t next_axb = row0[3] + row1[3]; - const int64_t sum_a = prev_a + next_a; - const int64_t sum_b = prev_b + next_b; - const int64_t sum_sq_ab = prev_sq + next_sq; - const int64_t sum_axb = prev_axb + next_axb; + const int32_t next_a = (int32_t)(row0[0] + row1[0]); + const int32_t next_b = (int32_t)(row0[1] + row1[1]); + const int32_t next_sq = (int32_t)(row0[2] + row1[2]); + const int32_t next_axb = (int32_t)(row0[3] + row1[3]); + const int32_t sum_a = prev_a + next_a; + const int32_t sum_b = prev_b + next_b; + const int32_t sum_sq_ab = prev_sq + next_sq; + const int32_t sum_axb = prev_axb + next_axb; prev_a = next_a; prev_b = next_b; prev_sq = next_sq; prev_axb = next_axb; - const int64_t sum_a_x_sum_b = sum_a * sum_b; - const int64_t ssim_n = (2 * sum_a_x_sum_b + c1) * - (2 * count * sum_axb - 2 * sum_a_x_sum_b + c2); - const int64_t sum_a_sq = sum_a * sum_a; - const int64_t sum_b_sq = sum_b * sum_b; - const int64_t ssim_d = (sum_a_sq + sum_b_sq + c1) * - (count * sum_sq_ab - sum_a_sq - sum_b_sq + c2); - ssim_total += ssim_n * 1.0 / ssim_d; + const int32_t sum_a_x_sum_b_2 = 2 * sum_a * sum_b; + const int32_t sum_a_sq_b_sq = sum_a * sum_a + sum_b * sum_b; + const int32_t n1 = sum_a_x_sum_b_2 + c1; + const int32_t n2 = (sum_axb << 7) - sum_a_x_sum_b_2 + c2; + const int32_t d1 = sum_a_sq_b_sq + c1; + const int32_t d2 = (sum_sq_ab << 6) - sum_a_sq_b_sq + c2; + ssim_total += ((float)n1 * (float)n2) / ((float)d1 * (float)d2); } return ssim_total; } diff --git a/source/compare_gcc.cc b/source/compare_gcc.cc index 2053e78c8..fcbb87de4 100644 --- a/source/compare_gcc.cc +++ b/source/compare_gcc.cc @@ -23,14 +23,20 @@ extern "C" { // asm statement, and other OS ABIs (e.g. Windows) differ in which registers // are callee-saved. List "memory" first, then "cc", then registers. -// This module is for GCC x86 and x64. -#if !defined(LIBYUV_DISABLE_X86) && \ - (defined(__x86_64__) || defined(__i386__)) && \ - !defined(LIBYUV_ENABLE_ROWWIN) +// This module is for GCC and clang x64. 32 bit x86 uses compare_win.cc. +#if !defined(LIBYUV_DISABLE_X86) && defined(__x86_64__) && \ + !defined(LIBYUV_ENABLE_ROWWIN) && \ + (!defined(_MSC_VER) || defined(__clang__)) + +#if defined(__clang__) || defined(__GNUC__) +#define LIBYUV_TARGET_AVX512BW \ + __attribute__((target("avx512bw,avx512vl,avx512f"))) +#else +#define LIBYUV_TARGET_AVX512BW +#endif // "memory" clobber prevents the reads from being removed -#if defined(__x86_64__) uint32_t HammingDistance_SSE42(const uint8_t* src_a, const uint8_t* src_b, int count) { @@ -78,46 +84,6 @@ uint32_t HammingDistance_SSE42(const uint8_t* src_a, return (uint32_t)(diff); } -#else -uint32_t HammingDistance_SSE42(const uint8_t* src_a, - const uint8_t* src_b, - int count) { - uint32_t diff = 0u; - - asm volatile( - // Process 16 bytes per loop. - LABELALIGN - "1: \n" - "mov (%0),%%ecx \n" - "mov 0x4(%0),%%edx \n" - "xor (%1),%%ecx \n" - "xor 0x4(%1),%%edx \n" - "popcnt %%ecx,%%ecx \n" - "add %%ecx,%3 \n" - "popcnt %%edx,%%edx \n" - "add %%edx,%3 \n" - "mov 0x8(%0),%%ecx \n" - "mov 0xc(%0),%%edx \n" - "xor 0x8(%1),%%ecx \n" - "xor 0xc(%1),%%edx \n" - "popcnt %%ecx,%%ecx \n" - "add %%ecx,%3 \n" - "popcnt %%edx,%%edx \n" - "add %%edx,%3 \n" - "add $0x10,%0 \n" - "add $0x10,%1 \n" - "sub $0x10,%2 \n" - "jg 1b \n" - : "+r"(src_a), // %0 - "+r"(src_b), // %1 - "+r"(count), // %2 - "+r"(diff) // %3 - : - : "cc", "memory", "ecx", "edx"); - - return diff; -} -#endif static const vec8 kNibbleMask = {15, 15, 15, 15, 15, 15, 15, 15, 15, 15, 15, 15, 15, 15, 15, 15}; @@ -419,6 +385,7 @@ static uint64_t LowBytesMask64(int n) { #ifdef HAS_HAMMINGDISTANCE_AVX512BW // Process 64 bytes per loop and the remainder of 1 to 63 bytes with masked // loads that read zeros, so any count is supported. +LIBYUV_TARGET_AVX512BW uint32_t HammingDistance_AVX512BW(const uint8_t* src_a, const uint8_t* src_b, int count) { @@ -487,6 +454,7 @@ uint32_t HammingDistance_AVX512BW(const uint8_t* src_a, #ifdef HAS_SUMSQUAREERROR_AVX512BW // Process 64 bytes per loop and the remainder of 1 to 63 bytes with masked // loads that read zeros, so any count is supported. +LIBYUV_TARGET_AVX512BW uint32_t SumSquareError_AVX512BW(const uint8_t* src_a, const uint8_t* src_b, int count) { @@ -704,6 +672,7 @@ static uint32_t Pow33(int n) { // Process 64 bytes per loop. The remainder of n = 1 to 63 bytes is loaded // into the high n bytes with a masked load ending at src + count, so the // same multipliers apply, and hash is multiplied by 33 ^ n. +LIBYUV_TARGET_AVX512BW uint32_t HashDjb2_AVX512BW(const uint8_t* src, int count, uint32_t seed) { const int n = count & 63; const uint64_t mask = n ? ~0ULL << (64 - n) : 0ULL; @@ -1043,6 +1012,7 @@ void Ssim4x4Row_AVX2(const uint8_t* src_a, // Process 64 pixels (16 4x4 tiles) per loop and the remainder of 4 to 60 // pixels (1 to 15 tiles) with masked loads and stores. +LIBYUV_TARGET_AVX512BW void Ssim4x4Row_AVX512BW(const uint8_t* src_a, int stride_a, const uint8_t* src_b, @@ -1121,7 +1091,257 @@ void Ssim4x4Row_AVX512BW(const uint8_t* src_a, #undef SSIM_TRANSPOSE_AVX512BW #endif // HAS_SSIM4X4ROW_AVX512BW -#endif // !defined(LIBYUV_DISABLE_X86) && (x86 or x64) && +#if defined(HAS_SSIMWINDOWROW_SSE2) || defined(HAS_SSIMWINDOWROW_AVX2) || \ + defined(HAS_SSIMWINDOWROW_AVX512BW) +static const uvec32 kSsimC1 = {26634, 26634, 26634, 26634}; +static const uvec32 kSsimC2 = {239708, 239708, 239708, 239708}; +#endif + +#ifdef HAS_SSIMWINDOWROW_SSE2 +// num_windows is a multiple of 4. +float SsimWindowRow_SSE2(const uint32_t* row0, + const uint32_t* row1, + int num_windows) { + float total; + asm volatile( + "sub %[row0],%[row1] \n" + "xorps %[total],%[total] \n" + + LABELALIGN + "1: \n" + "movdqu 0x00(%[row0]),%%xmm1 \n" + "movdqu 0x00(%[row0],%[row1]),%%xmm5 \n" + "paddd %%xmm5,%%xmm1 \n" + "movdqu 0x10(%[row0]),%%xmm2 \n" + "movdqu 0x10(%[row0],%[row1]),%%xmm5 \n" + "paddd %%xmm5,%%xmm2 \n" + "movdqu 0x20(%[row0]),%%xmm3 \n" + "movdqu 0x20(%[row0],%[row1]),%%xmm5 \n" + "paddd %%xmm5,%%xmm3 \n" + "movdqu 0x30(%[row0]),%%xmm4 \n" + "movdqu 0x30(%[row0],%[row1]),%%xmm5 \n" + "paddd %%xmm5,%%xmm4 \n" + "movdqu 0x40(%[row0]),%%xmm5 \n" + "add $0x40,%[row0] \n" + "paddd %%xmm2,%%xmm1 \n" + "paddd %%xmm3,%%xmm2 \n" + "paddd %%xmm4,%%xmm3 \n" + "paddd %%xmm5,%%xmm4 \n" + "movdqu 0x00(%[row0],%[row1]),%%xmm5 \n" + "paddd %%xmm5,%%xmm4 \n" + "movdqa %%xmm1,%%xmm5 \n" + "punpckldq %%xmm2,%%xmm1 \n" + "punpckhdq %%xmm2,%%xmm5 \n" + "movdqa %%xmm3,%%xmm2 \n" + "punpckldq %%xmm4,%%xmm3 \n" + "punpckhdq %%xmm4,%%xmm2 \n" + "movdqa %%xmm1,%%xmm4 \n" + "punpcklqdq %%xmm3,%%xmm1 \n" + "punpckhqdq %%xmm3,%%xmm4 \n" + "movdqa %%xmm5,%%xmm3 \n" + "punpcklqdq %%xmm2,%%xmm3 \n" + "punpckhqdq %%xmm2,%%xmm5 \n" + "movdqa %%xmm1,%%xmm2 \n" + "pmaddwd %%xmm4,%%xmm2 \n" + "pmaddwd %%xmm1,%%xmm1 \n" + "pmaddwd %%xmm4,%%xmm4 \n" + "paddd %%xmm4,%%xmm1 \n" + "pslld $0x1,%%xmm2 \n" + "pslld $0x7,%%xmm5 \n" + "psubd %%xmm2,%%xmm5 \n" + "paddd %[c1],%%xmm2 \n" + "paddd %[c2],%%xmm5 \n" + "pslld $0x6,%%xmm3 \n" + "psubd %%xmm1,%%xmm3 \n" + "paddd %[c1],%%xmm1 \n" + "paddd %[c2],%%xmm3 \n" + "cvtdq2ps %%xmm2,%%xmm2 \n" + "cvtdq2ps %%xmm5,%%xmm5 \n" + "cvtdq2ps %%xmm1,%%xmm1 \n" + "cvtdq2ps %%xmm3,%%xmm3 \n" + "mulps %%xmm5,%%xmm2 \n" + "mulps %%xmm3,%%xmm1 \n" + "divps %%xmm1,%%xmm2 \n" + "addps %%xmm2,%[total] \n" + "sub $0x4,%[num_windows] \n" + "jg 1b \n" + + "movhlps %[total],%%xmm1 \n" + "addps %%xmm1,%[total] \n" + "movaps %[total],%%xmm1 \n" + "shufps $0x1,%[total],%%xmm1 \n" + "addss %%xmm1,%[total] \n" + : [row0] "+r"(row0), [row1] "+r"(row1), [num_windows] "+r"(num_windows), + [total] "=&x"(total) + : [c1] "m"(kSsimC1), [c2] "m"(kSsimC2) + : "memory", "cc", "xmm1", "xmm2", "xmm3", "xmm4", "xmm5"); + return total; +} +#endif // HAS_SSIMWINDOWROW_SSE2 + +#ifdef HAS_SSIMWINDOWROW_AVX2 +// num_windows is a multiple of 8. +float SsimWindowRow_AVX2(const uint32_t* row0, + const uint32_t* row1, + int num_windows) { + float total; + asm volatile( + "sub %[row0],%[row1] \n" + "vxorps %t[total],%t[total],%t[total] \n" + + LABELALIGN + "1: \n" + "vmovdqu 0x00(%[row0]),%%ymm1 \n" + "vpaddd 0x00(%[row0],%[row1]),%%ymm1,%%ymm1 \n" + "vmovdqu 0x10(%[row0]),%%ymm5 \n" + "vpaddd 0x10(%[row0],%[row1]),%%ymm5,%%ymm5 \n" + "vpaddd %%ymm5,%%ymm1,%%ymm1 \n" + "vmovdqu 0x20(%[row0]),%%ymm2 \n" + "vpaddd 0x20(%[row0],%[row1]),%%ymm2,%%ymm2 \n" + "vmovdqu 0x30(%[row0]),%%ymm5 \n" + "vpaddd 0x30(%[row0],%[row1]),%%ymm5,%%ymm5 \n" + "vpaddd %%ymm5,%%ymm2,%%ymm2 \n" + "vmovdqu 0x40(%[row0]),%%ymm3 \n" + "vpaddd 0x40(%[row0],%[row1]),%%ymm3,%%ymm3 \n" + "vmovdqu 0x50(%[row0]),%%ymm5 \n" + "vpaddd 0x50(%[row0],%[row1]),%%ymm5,%%ymm5 \n" + "vpaddd %%ymm5,%%ymm3,%%ymm3 \n" + "vmovdqu 0x60(%[row0]),%%ymm4 \n" + "vpaddd 0x60(%[row0],%[row1]),%%ymm4,%%ymm4 \n" + "vmovdqu 0x70(%[row0]),%%ymm5 \n" + "vpaddd 0x70(%[row0],%[row1]),%%ymm5,%%ymm5 \n" + "vpaddd %%ymm5,%%ymm4,%%ymm4 \n" + "add $0x80,%[row0] \n" + "vpunpckldq %%ymm2,%%ymm1,%%ymm5 \n" + "vpunpckhdq %%ymm2,%%ymm1,%%ymm1 \n" + "vpunpckldq %%ymm4,%%ymm3,%%ymm2 \n" + "vpunpckhdq %%ymm4,%%ymm3,%%ymm3 \n" + "vpunpcklqdq %%ymm2,%%ymm5,%%ymm4 \n" + "vpunpckhqdq %%ymm2,%%ymm5,%%ymm2 \n" + "vpunpcklqdq %%ymm3,%%ymm1,%%ymm5 \n" + "vpunpckhqdq %%ymm3,%%ymm1,%%ymm3 \n" + "vpmaddwd %%ymm2,%%ymm4,%%ymm1 \n" + "vpmaddwd %%ymm4,%%ymm4,%%ymm4 \n" + "vpmaddwd %%ymm2,%%ymm2,%%ymm2 \n" + "vpaddd %%ymm2,%%ymm4,%%ymm4 \n" + "vpslld $0x1,%%ymm1,%%ymm1 \n" + "vpslld $0x7,%%ymm3,%%ymm3 \n" + "vpsubd %%ymm1,%%ymm3,%%ymm3 \n" + "vpslld $0x6,%%ymm5,%%ymm5 \n" + "vpsubd %%ymm4,%%ymm5,%%ymm5 \n" + "vpbroadcastd %[c1],%%ymm2 \n" + "vpaddd %%ymm2,%%ymm1,%%ymm1 \n" + "vpaddd %%ymm2,%%ymm4,%%ymm4 \n" + "vpbroadcastd %[c2],%%ymm2 \n" + "vpaddd %%ymm2,%%ymm3,%%ymm3 \n" + "vpaddd %%ymm2,%%ymm5,%%ymm5 \n" + "vcvtdq2ps %%ymm1,%%ymm1 \n" + "vcvtdq2ps %%ymm3,%%ymm3 \n" + "vcvtdq2ps %%ymm4,%%ymm4 \n" + "vcvtdq2ps %%ymm5,%%ymm5 \n" + "vmulps %%ymm3,%%ymm1,%%ymm1 \n" + "vmulps %%ymm5,%%ymm4,%%ymm4 \n" + "vdivps %%ymm4,%%ymm1,%%ymm1 \n" + "vaddps %%ymm1,%t[total],%t[total] \n" + "sub $0x8,%[num_windows] \n" + "jg 1b \n" + + "vextractf128 $0x1,%t[total],%%xmm1 \n" + "vaddps %%xmm1,%x[total],%x[total] \n" + "vmovhlps %x[total],%%xmm1,%%xmm1 \n" + "vaddps %%xmm1,%x[total],%x[total] \n" + "vshufps $0x1,%x[total],%x[total],%%xmm1 \n" + "vaddss %%xmm1,%x[total],%x[total] \n" + "vzeroupper \n" + : [row0] "+r"(row0), [row1] "+r"(row1), [num_windows] "+r"(num_windows), + [total] "=&x"(total) + : [c1] "m"(kSsimC1), [c2] "m"(kSsimC2) + : "memory", "cc", "xmm1", "xmm2", "xmm3", "xmm4", "xmm5"); + return total; +} +#endif // HAS_SSIMWINDOWROW_AVX2 + +#ifdef HAS_SSIMWINDOWROW_AVX512BW +// num_windows is a multiple of 16. +LIBYUV_TARGET_AVX512BW +float SsimWindowRow_AVX512BW(const uint32_t* row0, + const uint32_t* row1, + int num_windows) { + float total; + asm volatile( + "sub %[row0],%[row1] \n" + "vxorps %g[total],%g[total],%g[total] \n" + + LABELALIGN + "1: \n" + "vmovdqu32 0x00(%[row0]),%%zmm1 \n" + "vpaddd 0x00(%[row0],%[row1]),%%zmm1,%%zmm1 \n" + "vmovdqu32 0x40(%[row0]),%%zmm2 \n" + "vpaddd 0x40(%[row0],%[row1]),%%zmm2,%%zmm2 \n" + "vmovdqu32 0x80(%[row0]),%%zmm3 \n" + "vpaddd 0x80(%[row0],%[row1]),%%zmm3,%%zmm3 \n" + "vmovdqu32 0xc0(%[row0]),%%zmm4 \n" + "vpaddd 0xc0(%[row0],%[row1]),%%zmm4,%%zmm4 \n" + "valignd $0x4,%%zmm1,%%zmm2,%%zmm5 \n" + "vpaddd %%zmm5,%%zmm1,%%zmm1 \n" + "valignd $0x4,%%zmm2,%%zmm3,%%zmm5 \n" + "vpaddd %%zmm5,%%zmm2,%%zmm2 \n" + "valignd $0x4,%%zmm3,%%zmm4,%%zmm5 \n" + "vpaddd %%zmm5,%%zmm3,%%zmm3 \n" + "vmovdqu32 0xd0(%[row0]),%%zmm5 \n" + "vpaddd 0xd0(%[row0],%[row1]),%%zmm5,%%zmm5 \n" + "add $0x100,%[row0] \n" + "vpaddd %%zmm5,%%zmm4,%%zmm4 \n" + "vpunpckldq %%zmm2,%%zmm1,%%zmm5 \n" + "vpunpckhdq %%zmm2,%%zmm1,%%zmm1 \n" + "vpunpckldq %%zmm4,%%zmm3,%%zmm2 \n" + "vpunpckhdq %%zmm4,%%zmm3,%%zmm3 \n" + "vpunpcklqdq %%zmm2,%%zmm5,%%zmm4 \n" + "vpunpckhqdq %%zmm2,%%zmm5,%%zmm2 \n" + "vpunpcklqdq %%zmm3,%%zmm1,%%zmm5 \n" + "vpunpckhqdq %%zmm3,%%zmm1,%%zmm3 \n" + "vpmaddwd %%zmm2,%%zmm4,%%zmm1 \n" + "vpmaddwd %%zmm4,%%zmm4,%%zmm4 \n" + "vpmaddwd %%zmm2,%%zmm2,%%zmm2 \n" + "vpaddd %%zmm2,%%zmm4,%%zmm4 \n" + "vpslld $0x1,%%zmm1,%%zmm1 \n" + "vpslld $0x7,%%zmm3,%%zmm3 \n" + "vpsubd %%zmm1,%%zmm3,%%zmm3 \n" + "vpslld $0x6,%%zmm5,%%zmm5 \n" + "vpsubd %%zmm4,%%zmm5,%%zmm5 \n" + "vpaddd %[c1]%{1to16%},%%zmm1,%%zmm1 \n" + "vpaddd %[c2]%{1to16%},%%zmm3,%%zmm3 \n" + "vpaddd %[c1]%{1to16%},%%zmm4,%%zmm4 \n" + "vpaddd %[c2]%{1to16%},%%zmm5,%%zmm5 \n" + "vcvtdq2ps %%zmm1,%%zmm1 \n" + "vcvtdq2ps %%zmm3,%%zmm3 \n" + "vcvtdq2ps %%zmm4,%%zmm4 \n" + "vcvtdq2ps %%zmm5,%%zmm5 \n" + "vmulps %%zmm3,%%zmm1,%%zmm1 \n" + "vmulps %%zmm5,%%zmm4,%%zmm4 \n" + "vdivps %%zmm4,%%zmm1,%%zmm1 \n" + "vaddps %%zmm1,%g[total],%g[total] \n" + "sub $0x10,%[num_windows] \n" + "jg 1b \n" + + "vextractf32x8 $0x1,%g[total],%%ymm1 \n" + "vaddps %%ymm1,%t[total],%t[total] \n" + "vextractf128 $0x1,%t[total],%%xmm1 \n" + "vaddps %%xmm1,%x[total],%x[total] \n" + "vmovhlps %x[total],%%xmm1,%%xmm1 \n" + "vaddps %%xmm1,%x[total],%x[total] \n" + "vshufps $0x1,%x[total],%x[total],%%xmm1 \n" + "vaddss %%xmm1,%x[total],%x[total] \n" + "vzeroupper \n" + : [row0] "+r"(row0), [row1] "+r"(row1), [num_windows] "+r"(num_windows), + [total] "=&x"(total) + : [c1] "m"(kSsimC1), [c2] "m"(kSsimC2) + : "memory", "cc", "xmm1", "xmm2", "xmm3", "xmm4", "xmm5"); + return total; +} +#endif // HAS_SSIMWINDOWROW_AVX512BW + +#endif // !defined(LIBYUV_DISABLE_X86) && defined(__x86_64__) && // !defined(LIBYUV_ENABLE_ROWWIN) #ifdef __cplusplus diff --git a/source/compare_win.cc b/source/compare_win.cc index 469fa674d..c6a811901 100644 --- a/source/compare_win.cc +++ b/source/compare_win.cc @@ -11,13 +11,13 @@ #include "libyuv/compare_row.h" #include "libyuv/row.h" -// This module is for Visual C 32/64 bit, and for clang with -// LIBYUV_ENABLE_ROWWIN, which is also enabled for MemorySanitizer. -#if !defined(LIBYUV_DISABLE_X86) && \ - (defined(__x86_64__) || defined(__i386__) || defined(_M_X64) || \ - defined(_M_IX86)) && \ - ((defined(_MSC_VER) && !defined(__clang__)) || \ - defined(LIBYUV_ENABLE_ROWWIN)) +// This module is for 32 bit x86 with any compiler, Visual C 64 bit, and for +// clang with LIBYUV_ENABLE_ROWWIN, which is also enabled for MemorySanitizer. +#if !defined(LIBYUV_DISABLE_X86) && \ + (defined(__i386__) || defined(_M_IX86) || \ + ((defined(__x86_64__) || defined(_M_X64)) && \ + ((defined(_MSC_VER) && !defined(__clang__)) || \ + defined(LIBYUV_ENABLE_ROWWIN)))) #include // For SSE2 to AVX512 intrinsics #include // For memcpy @@ -37,7 +37,8 @@ extern "C" { #define LIBYUV_TARGET_SSE41 __attribute__((target("sse4.1"))) #define LIBYUV_TARGET_SSE42 __attribute__((target("sse4.2,popcnt"))) #define LIBYUV_TARGET_AVX2 __attribute__((target("avx2"))) -#define LIBYUV_TARGET_AVX512BW __attribute__((target("avx512bw,avx512f"))) +#define LIBYUV_TARGET_AVX512BW \ + __attribute__((target("avx512bw,avx512vl,avx512f"))) #else #define LIBYUV_TARGET_SSE2 #define LIBYUV_TARGET_SSE41 @@ -677,6 +678,194 @@ void Ssim4x4Row_AVX512BW(const uint8_t* src_a, } #endif // HAS_SSIM4X4ROW_AVX512BW +#if defined(HAS_SSIMWINDOWROW_SSE2) +// num_windows is a multiple of 4. +LIBYUV_TARGET_SSE2 +float SsimWindowRow_SSE2(const uint32_t* row0, + const uint32_t* row1, + int num_windows) { + const __m128i vc1 = _mm_set1_epi32(26634); + const __m128i vc2 = _mm_set1_epi32(239708); + __m128 total = _mm_setzero_ps(); + __m128i t0 = _mm_add_epi32(_mm_loadu_si128((const __m128i*)row0), + _mm_loadu_si128((const __m128i*)row1)); + int w; + for (w = 0; w < num_windows; w += 4) { + __m128i t1 = _mm_add_epi32(_mm_loadu_si128((const __m128i*)(row0 + 4)), + _mm_loadu_si128((const __m128i*)(row1 + 4))); + __m128i t2 = _mm_add_epi32(_mm_loadu_si128((const __m128i*)(row0 + 8)), + _mm_loadu_si128((const __m128i*)(row1 + 8))); + __m128i t3 = _mm_add_epi32(_mm_loadu_si128((const __m128i*)(row0 + 12)), + _mm_loadu_si128((const __m128i*)(row1 + 12))); + __m128i t4 = _mm_add_epi32(_mm_loadu_si128((const __m128i*)(row0 + 16)), + _mm_loadu_si128((const __m128i*)(row1 + 16))); + row0 += 16; + row1 += 16; + __m128i w0 = _mm_add_epi32(t0, t1); + __m128i w1 = _mm_add_epi32(t1, t2); + __m128i w2 = _mm_add_epi32(t2, t3); + __m128i w3 = _mm_add_epi32(t3, t4); + t0 = t4; + + __m128i u0 = _mm_unpacklo_epi32(w0, w1); + __m128i u1 = _mm_unpackhi_epi32(w0, w1); + __m128i u2 = _mm_unpacklo_epi32(w2, w3); + __m128i u3 = _mm_unpackhi_epi32(w2, w3); + __m128i sum_a = _mm_unpacklo_epi64(u0, u2); + __m128i sum_b = _mm_unpackhi_epi64(u0, u2); + __m128i sum_sq = _mm_unpacklo_epi64(u1, u3); + __m128i sum_axb = _mm_unpackhi_epi64(u1, u3); + + __m128i ab2 = _mm_slli_epi32(_mm_madd_epi16(sum_a, sum_b), 1); + __m128i a2b2 = _mm_add_epi32(_mm_madd_epi16(sum_a, sum_a), + _mm_madd_epi16(sum_b, sum_b)); + __m128i n1 = _mm_add_epi32(ab2, vc1); + __m128i n2 = + _mm_add_epi32(_mm_sub_epi32(_mm_slli_epi32(sum_axb, 7), ab2), vc2); + __m128i d1 = _mm_add_epi32(a2b2, vc1); + __m128i d2 = + _mm_add_epi32(_mm_sub_epi32(_mm_slli_epi32(sum_sq, 6), a2b2), vc2); + + __m128 num = _mm_mul_ps(_mm_cvtepi32_ps(n1), _mm_cvtepi32_ps(n2)); + __m128 den = _mm_mul_ps(_mm_cvtepi32_ps(d1), _mm_cvtepi32_ps(d2)); + total = _mm_add_ps(total, _mm_div_ps(num, den)); + } + total = _mm_add_ps(total, _mm_movehl_ps(total, total)); + total = _mm_add_ss(total, _mm_shuffle_ps(total, total, 1)); + return _mm_cvtss_f32(total); +} +#endif // HAS_SSIMWINDOWROW_SSE2 + +#if defined(HAS_SSIMWINDOWROW_AVX2) +// num_windows is a multiple of 8. +LIBYUV_TARGET_AVX2 +float SsimWindowRow_AVX2(const uint32_t* row0, + const uint32_t* row1, + int num_windows) { + const __m256i vc1 = _mm256_set1_epi32(26634); + const __m256i vc2 = _mm256_set1_epi32(239708); + __m256 total = _mm256_setzero_ps(); + __m128 total128; + int w; + for (w = 0; w < num_windows; w += 8) { + __m256i e01 = + _mm256_add_epi32(_mm256_loadu_si256((const __m256i*)(row0 + 0)), + _mm256_loadu_si256((const __m256i*)(row1 + 0))); + __m256i o12 = + _mm256_add_epi32(_mm256_loadu_si256((const __m256i*)(row0 + 4)), + _mm256_loadu_si256((const __m256i*)(row1 + 4))); + __m256i e23 = + _mm256_add_epi32(_mm256_loadu_si256((const __m256i*)(row0 + 8)), + _mm256_loadu_si256((const __m256i*)(row1 + 8))); + __m256i o34 = + _mm256_add_epi32(_mm256_loadu_si256((const __m256i*)(row0 + 12)), + _mm256_loadu_si256((const __m256i*)(row1 + 12))); + __m256i e45 = + _mm256_add_epi32(_mm256_loadu_si256((const __m256i*)(row0 + 16)), + _mm256_loadu_si256((const __m256i*)(row1 + 16))); + __m256i o56 = + _mm256_add_epi32(_mm256_loadu_si256((const __m256i*)(row0 + 20)), + _mm256_loadu_si256((const __m256i*)(row1 + 20))); + __m256i e67 = + _mm256_add_epi32(_mm256_loadu_si256((const __m256i*)(row0 + 24)), + _mm256_loadu_si256((const __m256i*)(row1 + 24))); + __m256i o78 = + _mm256_add_epi32(_mm256_loadu_si256((const __m256i*)(row0 + 28)), + _mm256_loadu_si256((const __m256i*)(row1 + 28))); + row0 += 32; + row1 += 32; + __m256i w0 = _mm256_add_epi32(e01, o12); + __m256i w1 = _mm256_add_epi32(e23, o34); + __m256i w2 = _mm256_add_epi32(e45, o56); + __m256i w3 = _mm256_add_epi32(e67, o78); + + __m256i u0 = _mm256_unpacklo_epi32(w0, w1); + __m256i u1 = _mm256_unpackhi_epi32(w0, w1); + __m256i u2 = _mm256_unpacklo_epi32(w2, w3); + __m256i u3 = _mm256_unpackhi_epi32(w2, w3); + __m256i sum_a = _mm256_unpacklo_epi64(u0, u2); + __m256i sum_b = _mm256_unpackhi_epi64(u0, u2); + __m256i sum_sq = _mm256_unpacklo_epi64(u1, u3); + __m256i sum_axb = _mm256_unpackhi_epi64(u1, u3); + + __m256i ab2 = _mm256_slli_epi32(_mm256_madd_epi16(sum_a, sum_b), 1); + __m256i a2b2 = _mm256_add_epi32(_mm256_madd_epi16(sum_a, sum_a), + _mm256_madd_epi16(sum_b, sum_b)); + __m256i n1 = _mm256_add_epi32(ab2, vc1); + __m256i n2 = _mm256_add_epi32( + _mm256_sub_epi32(_mm256_slli_epi32(sum_axb, 7), ab2), vc2); + __m256i d1 = _mm256_add_epi32(a2b2, vc1); + __m256i d2 = _mm256_add_epi32( + _mm256_sub_epi32(_mm256_slli_epi32(sum_sq, 6), a2b2), vc2); + + __m256 num = _mm256_mul_ps(_mm256_cvtepi32_ps(n1), _mm256_cvtepi32_ps(n2)); + __m256 den = _mm256_mul_ps(_mm256_cvtepi32_ps(d1), _mm256_cvtepi32_ps(d2)); + total = _mm256_add_ps(total, _mm256_div_ps(num, den)); + } + total128 = _mm_add_ps(_mm256_castps256_ps128(total), + _mm256_extractf128_ps(total, 1)); + total128 = _mm_add_ps(total128, _mm_movehl_ps(total128, total128)); + total128 = _mm_add_ss(total128, _mm_shuffle_ps(total128, total128, 1)); + _mm256_zeroupper(); + return _mm_cvtss_f32(total128); +} +#endif // HAS_SSIMWINDOWROW_AVX2 + +#if defined(HAS_SSIMWINDOWROW_AVX512BW) +// num_windows is a multiple of 16. +LIBYUV_TARGET_AVX512BW +float SsimWindowRow_AVX512BW(const uint32_t* row0, + const uint32_t* row1, + int num_windows) { + const __m512i vc1 = _mm512_set1_epi32(26634); + const __m512i vc2 = _mm512_set1_epi32(239708); + __m512 total = _mm512_setzero_ps(); + int w; + for (w = 0; w < num_windows; w += 16) { + __m512i t0 = _mm512_add_epi32(_mm512_loadu_si512(row0 + 0), + _mm512_loadu_si512(row1 + 0)); + __m512i t1 = _mm512_add_epi32(_mm512_loadu_si512(row0 + 16), + _mm512_loadu_si512(row1 + 16)); + __m512i t2 = _mm512_add_epi32(_mm512_loadu_si512(row0 + 32), + _mm512_loadu_si512(row1 + 32)); + __m512i t3 = _mm512_add_epi32(_mm512_loadu_si512(row0 + 48), + _mm512_loadu_si512(row1 + 48)); + __m512i t4 = _mm512_add_epi32(_mm512_loadu_si512(row0 + 52), + _mm512_loadu_si512(row1 + 52)); + row0 += 64; + row1 += 64; + __m512i w0 = _mm512_add_epi32(t0, _mm512_alignr_epi32(t1, t0, 4)); + __m512i w1 = _mm512_add_epi32(t1, _mm512_alignr_epi32(t2, t1, 4)); + __m512i w2 = _mm512_add_epi32(t2, _mm512_alignr_epi32(t3, t2, 4)); + __m512i w3 = _mm512_add_epi32(t3, t4); + + __m512i u0 = _mm512_unpacklo_epi32(w0, w1); + __m512i u1 = _mm512_unpackhi_epi32(w0, w1); + __m512i u2 = _mm512_unpacklo_epi32(w2, w3); + __m512i u3 = _mm512_unpackhi_epi32(w2, w3); + __m512i sum_a = _mm512_unpacklo_epi64(u0, u2); + __m512i sum_b = _mm512_unpackhi_epi64(u0, u2); + __m512i sum_sq = _mm512_unpacklo_epi64(u1, u3); + __m512i sum_axb = _mm512_unpackhi_epi64(u1, u3); + + __m512i ab2 = _mm512_slli_epi32(_mm512_madd_epi16(sum_a, sum_b), 1); + __m512i a2b2 = _mm512_add_epi32(_mm512_madd_epi16(sum_a, sum_a), + _mm512_madd_epi16(sum_b, sum_b)); + __m512i n1 = _mm512_add_epi32(ab2, vc1); + __m512i n2 = _mm512_add_epi32( + _mm512_sub_epi32(_mm512_slli_epi32(sum_axb, 7), ab2), vc2); + __m512i d1 = _mm512_add_epi32(a2b2, vc1); + __m512i d2 = _mm512_add_epi32( + _mm512_sub_epi32(_mm512_slli_epi32(sum_sq, 6), a2b2), vc2); + + __m512 num = _mm512_mul_ps(_mm512_cvtepi32_ps(n1), _mm512_cvtepi32_ps(n2)); + __m512 den = _mm512_mul_ps(_mm512_cvtepi32_ps(d1), _mm512_cvtepi32_ps(d2)); + total = _mm512_add_ps(total, _mm512_div_ps(num, den)); + } + return _mm512_reduce_add_ps(total); +} +#endif // HAS_SSIMWINDOWROW_AVX512BW + #ifdef __cplusplus } // extern "C" } // namespace libyuv diff --git a/unit_test/compare_test.cc b/unit_test/compare_test.cc index c6d19fb37..cfbe05144 100644 --- a/unit_test/compare_test.cc +++ b/unit_test/compare_test.cc @@ -664,16 +664,17 @@ TEST_F(LibYUVCompareTest, BenchmarkSsim_Opt) { MaskCpuFlags(benchmark_cpu_info_); - double opt_time = get_time(); + double ssim = 0.0; for (int i = 0; i < benchmark_iterations_; ++i) { - CalcFrameSsim(src_a, benchmark_width_, src_b, benchmark_width_, - benchmark_width_, benchmark_height_); + ssim = CalcFrameSsim(src_a, benchmark_width_, src_b, benchmark_width_, + benchmark_width_, benchmark_height_); } - opt_time = (get_time() - opt_time) / benchmark_iterations_; - printf("BenchmarkSsim_Opt - %8.2f us opt\n", opt_time * 1e6); - - ASSERT_EQ(0, 0); // Pass if we get this far. + if (benchmark_width_ > 8 && benchmark_height_ > 8) { + EXPECT_EQ(1.0, ssim); + } else { + EXPECT_EQ(0.0, ssim); + } free_aligned_buffer_page_end(src_a); free_aligned_buffer_page_end(src_b); @@ -689,16 +690,17 @@ TEST_F(LibYUVCompareTest, BenchmarkSsim_Unaligned) { MaskCpuFlags(benchmark_cpu_info_); - double opt_time = get_time(); + double ssim = 0.0; for (int i = 0; i < benchmark_iterations_; ++i) { - CalcFrameSsim(src_a + 1, benchmark_width_, src_b, benchmark_width_, - benchmark_width_, benchmark_height_); + ssim = CalcFrameSsim(src_a + 1, benchmark_width_, src_b, benchmark_width_, + benchmark_width_, benchmark_height_); } - opt_time = (get_time() - opt_time) / benchmark_iterations_; - printf("BenchmarkSsim_Opt - %8.2f us opt\n", opt_time * 1e6); - - ASSERT_EQ(0, 0); + if (benchmark_width_ > 8 && benchmark_height_ > 8) { + EXPECT_EQ(1.0, ssim); + } else { + EXPECT_EQ(0.0, ssim); + } free_aligned_buffer_page_end(src_a); free_aligned_buffer_page_end(src_b); @@ -710,8 +712,10 @@ TEST_F(LibYUVCompareTest, TestSsimAnyWidth) { const int kPlaneSize = kMaxWidth * kHeight + 1; align_buffer_page_end(src_a, kPlaneSize); align_buffer_page_end(src_b, kPlaneSize); - MemRandomize(src_a, kPlaneSize); - MemRandomize(src_b, kPlaneSize); + for (int i = 0; i < kPlaneSize; ++i) { + src_a[i] = (i * 73 + 19) & 0xff; + src_b[i] = (i * 151 + 43) & 0xff; + } for (int width = 9; width <= kMaxWidth; ++width) { for (int offset = 0; offset < 2; ++offset) { @@ -722,7 +726,8 @@ TEST_F(LibYUVCompareTest, TestSsimAnyWidth) { double c_ssim = CalcFrameSsim(a, width, b, width, width, kHeight); MaskCpuFlags(benchmark_cpu_info_); double opt_ssim = CalcFrameSsim(a, width, b, width, width, kHeight); - EXPECT_EQ(c_ssim, opt_ssim) << "width " << width << " offset " << offset; + EXPECT_NEAR(c_ssim, opt_ssim, 1e-5) + << "width " << width << " offset " << offset; } } @@ -745,15 +750,17 @@ TEST_F(LibYUVCompareTest, Ssim) { printf("warning - Ssim size too small. Testing function executes.\n"); } + // Identical images: SSIM == 1.0. double err; err = CalcFrameSsim(src_a + kSrcStride * b + b, kSrcStride, src_b + kSrcStride * b + b, kSrcStride, kSrcWidth, kSrcHeight); if (kSrcWidth > 8 && kSrcHeight > 8) { - ASSERT_EQ(err, 1.0); + EXPECT_EQ(err, 1.0); } + // Completely different constant images (255 vs 0): SSIM near 0. memset(src_a, 255, kSrcPlaneSize); err = CalcFrameSsim(src_a + kSrcStride * b + b, kSrcStride, @@ -761,9 +768,10 @@ TEST_F(LibYUVCompareTest, Ssim) { kSrcHeight); if (kSrcWidth > 8 && kSrcHeight > 8) { - ASSERT_LT(err, 0.0001); + EXPECT_LT(err, 0.0001); } + // Slightly different constant images (1 vs 0): SSIM ~ 0.84. memset(src_a, 1, kSrcPlaneSize); err = CalcFrameSsim(src_a + kSrcStride * b + b, kSrcStride, @@ -771,10 +779,11 @@ TEST_F(LibYUVCompareTest, Ssim) { kSrcHeight); if (kSrcWidth > 8 && kSrcHeight > 8) { - ASSERT_GT(err, 0.0001); - ASSERT_LT(err, 0.9); + EXPECT_GT(err, 0.8); + EXPECT_LT(err, 0.9); } + // Gradient vs constant 0: SSIM near 0. for (int i = 0; i < kSrcPlaneSize; ++i) { src_a[i] = i; } @@ -784,32 +793,37 @@ TEST_F(LibYUVCompareTest, Ssim) { kSrcHeight); if (kSrcWidth > 8 && kSrcHeight > 8) { - ASSERT_GT(err, 0.0); - ASSERT_LT(err, 0.01); + EXPECT_GT(err, 0.0); + EXPECT_LT(err, 0.01); } + // Similar images (deterministic gradient + small perturbation): SSIM ~ 0.99. for (int i = b; i < (kSrcHeight + b); ++i) { for (int j = b; j < (kSrcWidth + b); ++j) { - src_a[(i * kSrcStride) + j] = (fastrand() & 0xff); - src_b[(i * kSrcStride) + j] = (fastrand() & 0xff); + uint8_t val = static_cast((i * 3 + j * 5) & 0xff); + src_a[(i * kSrcStride) + j] = val; + src_b[(i * kSrcStride) + j] = + static_cast((val + ((i ^ j) & 3)) & 0xff); } } MaskCpuFlags(disable_cpu_flags_); - double c_err, opt_err; - - c_err = CalcFrameSsim(src_a + kSrcStride * b + b, kSrcStride, - src_b + kSrcStride * b + b, kSrcStride, kSrcWidth, - kSrcHeight); + double c_err = CalcFrameSsim(src_a + kSrcStride * b + b, kSrcStride, + src_b + kSrcStride * b + b, kSrcStride, + kSrcWidth, kSrcHeight); MaskCpuFlags(benchmark_cpu_info_); - - opt_err = CalcFrameSsim(src_a + kSrcStride * b + b, kSrcStride, - src_b + kSrcStride * b + b, kSrcStride, kSrcWidth, - kSrcHeight); + double opt_err = 0.0; + for (int i = 0; i < benchmark_iterations_; ++i) { + opt_err = CalcFrameSsim(src_a + kSrcStride * b + b, kSrcStride, + src_b + kSrcStride * b + b, kSrcStride, kSrcWidth, + kSrcHeight); + } if (kSrcWidth > 8 && kSrcHeight > 8) { - ASSERT_EQ(opt_err, c_err); + EXPECT_GT(opt_err, 0.9); + EXPECT_LT(opt_err, 1.0); + EXPECT_NEAR(c_err, opt_err, 1e-5); } free_aligned_buffer_page_end(src_a); diff --git a/util/compare.cc b/util/compare.cc index cd5b68dae..4ed933a1d 100644 --- a/util/compare.cc +++ b/util/compare.cc @@ -16,51 +16,126 @@ #include "libyuv/compare.h" #include "libyuv/version.h" +static bool ExtractResolutionFromFilename(const char* name, + int* width_ptr, + int* height_ptr) { + for (int i = 0; name[i]; ++i) { + if ((name[i] == '.' || name[i] == '_') && name[i + 1] >= '0' && + name[i + 1] <= '9') { + int n = sscanf(name + i + 1, "%dx%d", width_ptr, height_ptr); // NOLINT + if (2 == n) { + return true; + } + } + } + return false; +} + int main(int argc, char** argv) { - if (argc < 1) { + const char* name1 = NULL; + const char* name2 = NULL; + int width = 0; + int height = 0; + + for (int i = 1; i < argc; ++i) { + if (!strcmp(argv[i], "-s") && i + 2 < argc) { + width = atoi(argv[++i]); // NOLINT + height = atoi(argv[++i]); // NOLINT + } else if (!name1) { + name1 = argv[i]; + } else if (!name2) { + name2 = argv[i]; + } else if (width == 0 && i + 1 < argc) { + width = atoi(argv[i]); // NOLINT + height = atoi(argv[++i]); // NOLINT + } + } + + if (!name1) { printf("libyuv compare v%d\n", LIBYUV_VERSION); - printf("compare file1.yuv file2.yuv\n"); + printf("compare file1.yuv [file2.yuv [width height]]\n"); return -1; } - char* name1 = argv[1]; - char* name2 = (argc > 2) ? argv[2] : NULL; - FILE* fin1 = fopen(name1, "rb"); - FILE* fin2 = name2 ? fopen(name2, "rb") : NULL; - const int kBlockSize = 32768; - uint8_t buf1[kBlockSize]; - uint8_t buf2[kBlockSize]; + if (width <= 0 || height <= 0) { + if (!ExtractResolutionFromFilename(name1, &width, &height) && name2) { + ExtractResolutionFromFilename(name2, &width, &height); + } + } + + FILE* fin1 = fopen(name1, "rb"); + if (!fin1) { + fprintf(stderr, "Cannot open %s\n", name1); + return -1; + } + FILE* fin2 = name2 ? fopen(name2, "rb") : NULL; + if (name2 && !fin2) { + fprintf(stderr, "Cannot open %s\n", name2); + fclose(fin1); + return -1; + } + + const bool has_dims = (width >= 16 && height >= 16); + const int block_width = has_dims ? width : 256; + const int y_size = has_dims ? width * height : 32768; + const int uv_width = (width + 1) / 2; + const int uv_height = (height + 1) / 2; + const int uv_size = has_dims ? uv_width * uv_height : 0; + const int block_size = has_dims ? (y_size + 2 * uv_size) : 32768; + + uint8_t* buf1 = static_cast(malloc(block_size)); + uint8_t* buf2 = static_cast(malloc(block_size)); uint32_t hash1 = 5381; uint32_t hash2 = 5381; uint64_t sum_square_err = 0; uint64_t size_min = 0; + double ssim_sum = 0.0; + int ssim_count = 0; int amt1 = 0; int amt2 = 0; do { - amt1 = static_cast(fread(buf1, 1, kBlockSize, fin1)); + amt1 = static_cast(fread(buf1, 1, block_size, fin1)); if (amt1 > 0) { hash1 = libyuv::HashDjb2(buf1, amt1, hash1); } if (fin2) { - amt2 = static_cast(fread(buf2, 1, kBlockSize, fin2)); + amt2 = static_cast(fread(buf2, 1, block_size, fin2)); if (amt2 > 0) { hash2 = libyuv::HashDjb2(buf2, amt2, hash2); } int amt_min = (amt1 < amt2) ? amt1 : amt2; size_min += amt_min; sum_square_err += libyuv::ComputeSumSquareError(buf1, buf2, amt_min); + if (has_dims && amt_min == block_size) { + ssim_sum += libyuv::I420Ssim( + buf1, width, buf1 + y_size, uv_width, buf1 + y_size + uv_size, + uv_width, buf2, width, buf2 + y_size, uv_width, + buf2 + y_size + uv_size, uv_width, width, height); + ++ssim_count; + } else if (amt_min >= block_width * 8) { + ssim_sum += libyuv::CalcFrameSsim(buf1, block_width, buf2, block_width, + block_width, amt_min / block_width); + ++ssim_count; + } } } while (amt1 > 0 || amt2 > 0); printf("hash1 %x", hash1); if (fin2) { printf(", hash2 %x", hash2); - double mse = - static_cast(sum_square_err) / static_cast(size_min); + double mse = (size_min > 0) ? static_cast(sum_square_err) / + static_cast(size_min) + : 0.0; printf(", mse %.2f", mse); double psnr = libyuv::SumSquareErrorToPsnr(sum_square_err, size_min); - printf(", psnr %.2f\n", psnr); + printf(", psnr %.2f", psnr); + double ssim = (ssim_count > 0) ? (ssim_sum / ssim_count) : 0.0; + printf(", ssim %.6f", ssim); fclose(fin2); } + printf("\n"); fclose(fin1); + free(buf1); + free(buf2); + return 0; }