Commit ee674323 authored by David Barker's avatar David Barker Committed by Debargha Mukherjee

Vectorize corner matching function

Add an SSE4 version of compute_cross_correlation() from
corner_match.c. This function is about 3.4x the speed of
the scalar code; determine_correspondence as a whole is about
2.5-3x the speed it was previously.

BUG=aomedia:487

Change-Id: I707b7cfd5c513c025d3ee7fb6a5f1fa335ecd495
parent f695d659
......@@ -163,4 +163,8 @@ AV1_CX_SRCS-$(HAVE_MSA) += encoder/mips/msa/fdct16x16_msa.c
AV1_CX_SRCS-$(HAVE_MSA) += encoder/mips/msa/fdct_msa.h
AV1_CX_SRCS-$(HAVE_MSA) += encoder/mips/msa/temporal_filter_msa.c
ifeq ($(CONFIG_GLOBAL_MOTION),yes)
AV1_CX_SRCS-$(HAVE_SSE4_1) += encoder/x86/corner_match_sse4.c
endif
AV1_CX_SRCS-yes := $(filter-out $(AV1_CX_SRCS_REMOVE-yes),$(AV1_CX_SRCS-yes))
......@@ -617,6 +617,12 @@ if ((aom_config("CONFIG_WARPED_MOTION") eq "yes") ||
}
}
if (aom_config("CONFIG_WARPED_MOTION") eq "yes" &&
aom_config("CONFIG_AV1_ENCODER") eq "yes") {
add_proto qw/double compute_cross_correlation/, "unsigned char *im1, int stride1, int x1, int y1, unsigned char *im2, int stride2, int x2, int y2";
specialize qw/compute_cross_correlation sse4_1/;
}
# LOOP_RESTORATION functions
if (aom_config("CONFIG_LOOP_RESTORATION") eq "yes") {
......
......@@ -9,16 +9,13 @@
* PATENTS file, you can obtain it at www.aomedia.org/license/patent.
*/
#include <stdio.h>
#include <stdlib.h>
#include <memory.h>
#include <math.h>
#include "./av1_rtcd.h"
#include "av1/encoder/corner_match.h"
#define MATCH_SZ 13
#define MATCH_SZ_BY2 ((MATCH_SZ - 1) / 2)
#define MATCH_SZ_SQ (MATCH_SZ * MATCH_SZ)
#define SEARCH_SZ 9
#define SEARCH_SZ_BY2 ((SEARCH_SZ - 1) / 2)
......@@ -46,9 +43,9 @@ static double compute_variance(unsigned char *im, int stride, int x, int y) {
correlation/standard deviation are taken over MATCH_SZ by MATCH_SZ windows
of each image, centered at (x1, y1) and (x2, y2) respectively.
*/
static double compute_cross_correlation(unsigned char *im1, int stride1, int x1,
int y1, unsigned char *im2, int stride2,
int x2, int y2) {
double compute_cross_correlation_c(unsigned char *im1, int stride1, int x1,
int y1, unsigned char *im2, int stride2,
int x2, int y2) {
int v1, v2;
int sum1 = 0;
int sum2 = 0;
......
......@@ -15,6 +15,10 @@
#include <stdlib.h>
#include <memory.h>
#define MATCH_SZ 13
#define MATCH_SZ_BY2 ((MATCH_SZ - 1) / 2)
#define MATCH_SZ_SQ (MATCH_SZ * MATCH_SZ)
typedef struct {
int x, y;
int rx, ry;
......
#include <stdlib.h>
#include <memory.h>
#include <math.h>
#include <assert.h>
#include <smmintrin.h>
#include "./av1_rtcd.h"
#include "aom_ports/mem.h"
#include "av1/encoder/corner_match.h"
DECLARE_ALIGNED(16, static const uint8_t, byte_mask[16]) = {
255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 0, 0, 0
};
#if MATCH_SZ != 13
#error "Need to change byte_mask in corner_match_sse4.c if MATCH_SZ != 13"
#endif
/* Compute corr(im1, im2) * MATCH_SZ * stddev(im1), where the
correlation/standard deviation are taken over MATCH_SZ by MATCH_SZ windows
of each image, centered at (x1, y1) and (x2, y2) respectively.
*/
double compute_cross_correlation_sse4_1(unsigned char *im1, int stride1, int x1,
int y1, unsigned char *im2, int stride2,
int x2, int y2) {
int i;
// 2 16-bit partial sums in lanes 0, 4 (== 2 32-bit partial sums in lanes 0,
// 2)
__m128i sum1_vec = _mm_setzero_si128();
__m128i sum2_vec = _mm_setzero_si128();
// 4 32-bit partial sums of squares
__m128i sumsq2_vec = _mm_setzero_si128();
__m128i cross_vec = _mm_setzero_si128();
const __m128i mask = _mm_load_si128((__m128i *)byte_mask);
const __m128i zero = _mm_setzero_si128();
im1 += (y1 - MATCH_SZ_BY2) * stride1 + (x1 - MATCH_SZ_BY2);
im2 += (y2 - MATCH_SZ_BY2) * stride2 + (x2 - MATCH_SZ_BY2);
for (i = 0; i < MATCH_SZ; ++i) {
const __m128i v1 =
_mm_and_si128(_mm_loadu_si128((__m128i *)&im1[i * stride1]), mask);
const __m128i v2 =
_mm_and_si128(_mm_loadu_si128((__m128i *)&im2[i * stride2]), mask);
// Using the 'sad' intrinsic here is a bit faster than adding
// v1_l + v1_r and v2_l + v2_r, plus it avoids the need for a 16->32 bit
// conversion step later, for a net speedup of ~10%
sum1_vec = _mm_add_epi16(sum1_vec, _mm_sad_epu8(v1, zero));
sum2_vec = _mm_add_epi16(sum2_vec, _mm_sad_epu8(v2, zero));
const __m128i v1_l = _mm_cvtepu8_epi16(v1);
const __m128i v1_r = _mm_cvtepu8_epi16(_mm_srli_si128(v1, 8));
const __m128i v2_l = _mm_cvtepu8_epi16(v2);
const __m128i v2_r = _mm_cvtepu8_epi16(_mm_srli_si128(v2, 8));
sumsq2_vec = _mm_add_epi32(
sumsq2_vec,
_mm_add_epi32(_mm_madd_epi16(v2_l, v2_l), _mm_madd_epi16(v2_r, v2_r)));
cross_vec = _mm_add_epi32(
cross_vec,
_mm_add_epi32(_mm_madd_epi16(v1_l, v2_l), _mm_madd_epi16(v1_r, v2_r)));
}
// Now we can treat the four registers (sum1_vec, sum2_vec, sumsq2_vec,
// cross_vec)
// as holding 4 32-bit elements each, which we want to sum horizontally.
// We do this by transposing and then summing vertically.
__m128i tmp_0 = _mm_unpacklo_epi32(sum1_vec, sum2_vec);
__m128i tmp_1 = _mm_unpackhi_epi32(sum1_vec, sum2_vec);
__m128i tmp_2 = _mm_unpacklo_epi32(sumsq2_vec, cross_vec);
__m128i tmp_3 = _mm_unpackhi_epi32(sumsq2_vec, cross_vec);
__m128i tmp_4 = _mm_unpacklo_epi64(tmp_0, tmp_2);
__m128i tmp_5 = _mm_unpackhi_epi64(tmp_0, tmp_2);
__m128i tmp_6 = _mm_unpacklo_epi64(tmp_1, tmp_3);
__m128i tmp_7 = _mm_unpackhi_epi64(tmp_1, tmp_3);
__m128i res =
_mm_add_epi32(_mm_add_epi32(tmp_4, tmp_5), _mm_add_epi32(tmp_6, tmp_7));
int sum1 = _mm_extract_epi32(res, 0);
int sum2 = _mm_extract_epi32(res, 1);
int sumsq2 = _mm_extract_epi32(res, 2);
int cross = _mm_extract_epi32(res, 3);
int var2 = sumsq2 * MATCH_SZ_SQ - sum2 * sum2;
int cov = cross * MATCH_SZ_SQ - sum1 * sum2;
return cov / sqrt((double)var2);
}
/*
* Copyright (c) 2016, Alliance for Open Media. All rights reserved
*
* This source code is subject to the terms of the BSD 2 Clause License and
* the Alliance for Open Media Patent License 1.0. If the BSD 2 Clause License
* was not distributed with this source code in the LICENSE file, you can
* obtain it at www.aomedia.org/license/software. If the Alliance for Open
* Media Patent License 1.0 was not distributed with this source code in the
* PATENTS file, you can obtain it at www.aomedia.org/license/patent.
*/
#include "third_party/googletest/src/googletest/include/gtest/gtest.h"
#include "test/acm_random.h"
#include "test/util.h"
#include "./av1_rtcd.h"
#include "test/clear_system_state.h"
#include "test/register_state_check.h"
#include "av1/encoder/corner_match.h"
namespace test_libaom {
namespace AV1CornerMatch {
using libaom_test::ACMRandom;
using std::tr1::tuple;
using std::tr1::make_tuple;
typedef tuple<int> CornerMatchParam;
class AV1CornerMatchTest : public ::testing::TestWithParam<CornerMatchParam> {
public:
virtual ~AV1CornerMatchTest();
virtual void SetUp();
virtual void TearDown();
protected:
void RunCheckOutput();
libaom_test::ACMRandom rnd_;
};
AV1CornerMatchTest::~AV1CornerMatchTest() {}
void AV1CornerMatchTest::SetUp() { rnd_.Reset(ACMRandom::DeterministicSeed()); }
void AV1CornerMatchTest::TearDown() { libaom_test::ClearSystemState(); }
void AV1CornerMatchTest::RunCheckOutput() {
const int w = 128, h = 128;
const int num_iters = 10000;
int i, j;
uint8_t *input1 = new uint8_t[w * h];
uint8_t *input2 = new uint8_t[w * h];
// Test the two extreme cases:
// i) Random data, should have correlation close to 0
// ii) Linearly related data + noise, should have correlation close to 1
int mode = GET_PARAM(0);
if (mode == 0) {
for (i = 0; i < h; ++i)
for (j = 0; j < w; ++j) {
input1[i * w + j] = rnd_.Rand8();
input2[i * w + j] = rnd_.Rand8();
}
} else if (mode == 1) {
for (i = 0; i < h; ++i)
for (j = 0; j < w; ++j) {
int v = rnd_.Rand8();
input1[i * w + j] = v;
input2[i * w + j] = (v / 2) + (rnd_.Rand8() & 15);
}
}
for (i = 0; i < num_iters; ++i) {
int x1 = MATCH_SZ_BY2 + rnd_.PseudoUniform(w - 2 * MATCH_SZ_BY2);
int y1 = MATCH_SZ_BY2 + rnd_.PseudoUniform(h - 2 * MATCH_SZ_BY2);
int x2 = MATCH_SZ_BY2 + rnd_.PseudoUniform(w - 2 * MATCH_SZ_BY2);
int y2 = MATCH_SZ_BY2 + rnd_.PseudoUniform(h - 2 * MATCH_SZ_BY2);
double res_c =
compute_cross_correlation_c(input1, w, x1, y1, input2, w, x2, y2);
double res_sse4 =
compute_cross_correlation_sse4_1(input1, w, x1, y1, input2, w, x2, y2);
ASSERT_EQ(res_sse4, res_c);
}
delete[] input1;
delete[] input2;
}
TEST_P(AV1CornerMatchTest, CheckOutput) { RunCheckOutput(); }
INSTANTIATE_TEST_CASE_P(SSE4_1, AV1CornerMatchTest,
::testing::Values(make_tuple(0), make_tuple(1)));
} // namespace AV1CornerMatch
} // namespace test_libaom
......@@ -233,6 +233,10 @@ ifeq ($(CONFIG_LOOP_RESTORATION),yes)
LIBAOM_TEST_SRCS-$(HAVE_SSE4_1) += selfguided_filter_test.cc
endif
ifeq ($(CONFIG_GLOBAL_MOTION)$(CONFIG_AV1_ENCODER),yesyes)
LIBAOM_TEST_SRCS-$(HAVE_SSE4_1) += corner_match_test.cc
endif
TEST_INTRA_PRED_SPEED_SRCS-yes := test_intra_pred_speed.cc
TEST_INTRA_PRED_SPEED_SRCS-yes += ../md5_utils.h ../md5_utils.c
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment