优化 2D 旋转

Question

给定二维旋转点的经典公式space：

cv::Point pt[NPOINTS];
cv::Point rotated[NPOINTS];
float angle = WHATEVER;
float cosine = cos(angle);
float sine = sin(angle);

for (int i = 0; i < NPOINTS; i++)
{
    rotated[i].x = pt[i].x * cosine - pt[i].y * sine;
    rotated[i].y = pt[i].x * sine   + pt[i].y * cosine;
}

给定 NPOINTS 为 32 且数组对齐，如何优化 SSE 或 AVX 的代码？在这里和其他地方搜索没有找到任何有用的东西，我在这里迷路了：

__m128i onePoint = _mm_set_epi32(pt[i].x, pt[i].y, pt[i].x, pt[i].y);
__m128 onefPoint = _m128_cvtepi32_ps(onePoint);
__m128 sinCos = _mm_set_ps(cosine, -sine, sine, cosine);
__m128 rotated = _mm_mul_ps(onefPoint, sinCos);

但是如何从[y*cosine, -x*sine, x*sine, y*cosine]到[y*cosine + -x*sine, x*sine + y*cosine]呢？这是最好的方法吗？它可以轻松扩展到 __m512 吗？

更新：我做了更多的研究，现在大约有：

__m128i onePoint = _mm_set_epi32(pt[i].x, pt[i].y, pt[i].x, pt[i].y);
__m128 onefPoint = _m128_cvtepi32_ps(onePoint);
__m128i twoPoint = _mm_set_epi32(pt[i+1].x, pt[i+1].y, pt[i+1].x, pt[i+1].y);
__m128 twofPoint = _m128_cvtepi32_ps(twoPoint);
__m128 sinCos = _mm_set_ps(cosine, -sine, sine, cosine);
__m128 rotated1 = _mm_mul_ps(onefPoint, sinCos);
__m128 rotated2 = _mm_mul_ps(twofPoint, sinCos);
__m128 added = _mm_hadd_ps(rotated1, rotated2);
__m128i intResult = _mm_cvtps_epi32(added);
int results[4];
_mm_storeu_si128((__m128i*)results, intResult);

这将 50% 的速度从 11% 的处理器时间提高到大约 6%。扩展到 __m256 并一次做四个点可以提供另一个加速。这看起来很糟糕的代码，但我的方向正确吗？

Answer 1

使用数组结构数组 (AoSoA) 并一次处理八个点。在下面的代码中，point8 是包含八个点的数组结构。函数 rotate_point8 旋转八个点，与旋转一个点的函数 rotate_point 具有相同的代数结构。函数 rotate_all8 使用 AoSoA point8*.

旋转 32 个点

单点旋转代码做了4次乘法1次加法1次减法

如果我们查看 the assembly for rotate_point8，我们会看到 GCC 展开循环并在每次展开时执行 4 次 SIMD 乘法、1 次 SIMD 加法和 1 次 SIMD 减法。那是你能做的最好的事情：八个一个的价格。

#include <x86intrin.h>
#include <stdio.h>
#include <math.h>

struct point8 {
  __m256 x;
  __m256 y;
};

struct point {
  float x;
  float y;
};

static point rotate_point(point p, float a, float b) {
  point r;
  r.x = p.x*a - p.y*b;
  r.y = p.x*b + p.y*a;
  return r;
}

static point8 rotate_point8(point8 p, float a, float b) {
  __m256 va = _mm256_set1_ps(a), vb = _mm256_set1_ps(b);
  point8 r;
  r.x = _mm256_sub_ps(_mm256_mul_ps(p.x,va), _mm256_mul_ps(p.y,vb));
  r.y = _mm256_add_ps(_mm256_mul_ps(p.x,vb), _mm256_mul_ps(p.y,va));
  return r;
}

void rotate_all(point* points, point* r, float angle) {
  float a = cos(angle), b = sin(angle);
  for(int i=0; i<32; i++) r[i] = rotate_point(points[i], a, b);
}

void rotate_all8(point8* points, point8* r8, float angle) {
  float a = cos(angle), b = sin(angle);
  for(int i=0; i<4; i++) r8[i] = rotate_point8(points[i], a, b);
}

int main(void) {
  float x[32], y[32];
  point p[32], r[32];
  point8 p8[4], r8[4];
  float angle = 3.14159f/4;

  for(int i=0; i<32; i++) y[i] = 1.0*i/31, x[i] = sqrt(1-y[i]*y[i]);
  for(int i=0; i<32; i++) p[i].x = x[i], p[i].y = y[i];
  for(int i=0; i<4; i++) p8[i].x = _mm256_load_ps(&x[8*i]), p8[i].y = _mm256_load_ps(&y[8*i]); 

  for(int i=0; i<32; i++) printf("%f %f\n", p[i].x, p[i].y); puts("");

  rotate_all(p, r, angle);
  for(int i=0; i<32; i++) printf("%f %f\n", r[i].x, r[i].y); puts("");

  rotate_all8(p8, r8, angle);
  for(int i=0; i<4; i++) {
    _mm256_storeu_ps(x, r8[i].x), _mm256_storeu_ps(y, r8[i].y);
    for(int j=0; j<8; j++) printf("%f %f\n", x[j], y[j]);
  }
}

优化 2D 旋转

Optimising 2D rotation

c++

optimization

opencv

avx