检测点是否在扇形之内

最新推荐文章于 2024-01-01 09:35:06 发布
最新推荐文章于 2024-01-01 09:35:06 发布
阅读量2.1k 收藏 1
点赞数
分类专栏： 图形算法与图形硬件

前几天，同事在报告中提及检测角色是否在扇形攻击范围的方法。我觉得该方法的性能不是太好，提出另一个颇为直接的方法。此问题在游戏中十分常见，只涉及简单的数学，却又可以看出实现者是否细心，所以我觉得可当作一道简单的面试题。问题在微博发表后得到不少回应，故撰文提供一些解答。

问题定义:

在二维中，检测点 是否在 扇形(circular sector)内，设扇形的顶点为 ，半径为 ，从 方向两边展开 角度。

当中  以直角坐标(cartesian coordinates)表示，，。在扇形区边界上当作不相交。实现时所有数值采用单精度浮点数类型。

问题分析

许多相交问题都可以分拆为较小的问题。在此问题中，扇形可表示为圆形和角度区间的交集。

换句话说，问题等同于检测  和  的距离小于 ，及  的方向在  两边  的角度范围内。

距离

 和  的距离小于 ， 用数学方式表示:

极坐标

这是比较麻烦的部分，也是本题的核心。

有人想到，可以把  和  从直角坐标转换成极坐标(polar coordinates)。数学上， 和  分别与  轴的夹角可用atan2()函数求得:

然后，检查  是否在  区间内。但这要非常小心，因为  区间可能超越  的范围，所以要检测:

这个方法是可行的，不过即使假设  和  是常数，可预计算  和  ，我们还是避免不了要计算一个atan2()。

点积

点积(dot product)可计算两个矢量的夹角，这非常适合本题的扇形表示方式。我们要检测  和  的夹角是否小于  :

相比极坐标的方法，点积算出来的夹角必然在  区间里，无需作特别处理就可以和  比较。

这是比较直观的角度比较方式。本文将以此方法为主轴。

编码与优化

若直接实现以距离和点积的检测，可以得到:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
// Naive
bool  IsPointInCircularSector(
     float  cx, float  cy, float  ux, float  uy, float  r, float  theta,
     float  px, float  py)
{
     assert (cosTheta > -1 && cosTheta < 1);
     assert (squaredR > 0.0f);
 
     // D = P - C
     float  dx = px - cx;
     float  dy = py - cy;
 
     // |D| = (dx^2 + dy^2)^0.5
     float  length = sqrt (dx * dx + dy * dy);
 
     // |D| > r
     if  (length > r)
         return  false ;
 
     // Normalize D
     dx /= length;
     dy /= length;
 
     // acos(D dot U) < theta
     return  acos (dx * ux + dy * uy) < theta;
}

优化版本1: 基本优化

由于计算矢量长度需要sqrt()，而不等式(1.1)左右两侧都是非负数，所以可以把两侧平方:

如果  为常数，我们可以预计算  。

另外，如果是常数，我们可以预计算  ，然后把不等式(1.5) 改为:

可以这么做，是因为  在  的区间内是单调下降的。

此外，由于除法一般较乘法慢得多，而  ，我们可以把  调至右侧:

基于这两个可能预计算的参数，可把检测函数的参数由  和  改成  和  。结合这些改动:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
// Basic: use squareR and cosTheta as parameters, defer sqrt(), eliminate division
bool  IsPointInCircularSector1(
     float  cx, float  cy, float  ux, float  uy, float  squaredR, float  cosTheta,
     float  px, float  py)
{
     assert (cosTheta > -1 && cosTheta < 1);
     assert (squaredR > 0.0f);
 
     // D = P - C
     float  dx = px - cx;
     float  dy = py - cy;
 
     // |D|^2 = (dx^2 + dy^2)
     float  squaredLength = dx * dx + dy * dy;
 
     // |D|^2 > r^2
     if  (squaredLength > squaredR)
         return  false ;
 
     // |D|
     float  length = sqrt (squaredLength);
 
     // D dot U > |D| cos(theta)
     return  dx * ux + dy * uy > length * cosTheta;
}

注意，虽然比较长度时不用开平方，在夹角的检测里还是要算一次开平方，但这也比必须算开平方好，因为第一个检测失败就不用算了。

优化版本2: 去除开平方

然而，我们有办法去除开平方吗?

这或许是这解答中最难的问题。我们不能简单地把不等式(2.3)两侧平方:

因为两侧分别有可能是正数或负数。若把负数的一侧平方后，就会把负号清除了，我们必须把比较符号反转。

针对左右侧分别可以是负数和非负数四种情况，可以归纳为:

  1. 若左侧为非负数，右侧为非负数，检测 
  2. 若左侧为为负数，右侧为负数，检测 
  3. 若左侧为非负数，右侧为负数，那么不等式(2.3)一定成立
  4. 若左侧为负数，右侧为非负数，那么不等式(2.3)一定不成立

解释一下，在第2个情况中，先把两侧分别乘以-1，大于就变成小于，两侧就变成非负，可以平方。在实现中，若非第1个或第2个情况，只可能是第3个或第4个，所以只需分辨是3还是4，例如只检测左侧是否负数。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
// Eliminate sqrt()
bool  IsPointInCircularSector2(
     float  cx, float  cy, float  ux, float  uy, float  squaredR, float  cosTheta,
     float  px, float  py)
{
     assert (cosTheta > -1 && cosTheta < 1);
     assert (squaredR > 0.0f);
 
     // D = P - C
     float  dx = px - cx;
     float  dy = py - cy;
 
     // |D|^2 = (dx^2 + dy^2)
     float  squaredLength = dx * dx + dy * dy;
 
     // |D|^2 > r^2
     if  (squaredLength > squaredR)
         return  false ;
 
     // D dot U
     float  DdotU = dx * ux + dy * uy;
 
     // D dot U > |D| cos(theta)
     // <=>
     // (D dot U)^2 > |D|^2 (cos(theta))^2 if D dot U >= 0 and cos(theta) >= 0
     // (D dot U)^2 < |D|^2 (cos(theta))^2 if D dot U <  0 and cos(theta) <  0
     // true                               if D dot U >= 0 and cos(theta) <  0
     // false                              if D dot U <  0 and cos(theta) >= 0
     if  (DdotU >= 0 && cosTheta >= 0)
         return  DdotU * DdotU > squaredLength * cosTheta * cosTheta;
     else  if  (DdotU < 0 && cosTheta < 0)
         return  DdotU * DdotU < squaredLength * cosTheta * cosTheta;
     else
         return  DdotU >= 0;
}

优化版本3: 减少比较

在一些架构上，浮点比较和分支是较慢的操作。由于我们只是要检测数值是否负值，我们可以利用IEEE754浮点数的二进制特性：最高位代表符号，0为正数，1为负数。

我们可把符号用位掩码(bit mask)分离出来，第1和第2个情况就变成比较左右侧的符号是否相等。我们还可以用OR把符号加到两侧平方，意义上就是当负数时把比较符号调换，

而第3个情况也可以直接用符号作比较。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
// Bit trick
bool  IsPointInCircularSector3(
     float  cx, float  cy, float  ux, float  uy, float  squaredR, float  cosTheta,
     float  px, float  py)
{
     assert (cosTheta > -1 && cosTheta < 1);
     assert (squaredR > 0.0f);
 
     // D = P - C
     float  dx = px - cx;
     float  dy = py - cy;
 
     // |D|^2 = (dx^2 + dy^2)
     float  squaredLength = dx * dx + dy * dy;
 
     // |D|^2 > r^2
     if  (squaredLength > squaredR)
         return  false ;
 
     // D dot U
     float  DdotU = dx * ux + dy * uy;
 
     // D dot U > |D| cos(theta)
     // <=>
     // (D dot U)^2 > |D|^2 (cos(theta))^2 if D dot U >= 0 and cos(theta) >= 0
     // (D dot U)^2 < |D|^2 (cos(theta))^2 if D dot U <  0 and cos(theta) <  0
     // true                               if D dot U >= 0 and cos(theta) <  0
     // false                              if D dot U <  0 and cos(theta) >= 0
     const  unsigned cSignMask = 0x80000000;
     union  {
         float  f;
         unsigned u;
     }a, b, lhs, rhs;
     a.f = DdotU;
     b.f = cosTheta;
     unsigned asign = a.u & cSignMask;
     unsigned bsign = b.u & cSignMask;
     if  (asign == bsign) {
         lhs.f = DdotU * DdotU;
         rhs.f = squaredLength * cosTheta * cosTheta;
         lhs.u |= asign;
         rhs.u |= asign;
         return  lhs.f > rhs.f;
     }
     else
         return  asign == 0;
}
<h2<性能测试< h2="" style="color: rgb(51, 51, 51); font-family: Georgia, 'Times New Roman', Times, sans-serif; font-size: 14px; line-height: 25.2px;">

上述所作的”优化”，其实只是从算法上著手，试图消除一些可能开销较高的操作。实际上不同的CPU架构的结果也会有所差异。所以我们必须实测才能得知，在某硬件、编译器上，”优化”是否有效。

时间所限，我只做了一个简单的性能测试(欠缺足够边界条件单元测试……):

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
const  unsigned N = 1000;
const  unsigned M = 100000;
 
template  < bool  naiveParam, typename  TestFunc>
float  Test(TestFunc f, float  rmax = 2.0f) {
     unsigned count = 0;
     for  ( int  i = 0; i < N; i++) {
         float  cx = gCx[i];
         float  cy = gCy[i];
         float  ux = gUx[i];
         float  uy = gUy[i];
         float  r = naiveParam ? gR[i] : gSquaredR[i];
         float  t = naiveParam ? gTheta[i] : gCosTheta[i];
 
         for  ( int  j = 0; j < M; j++) {
             if  (f(cx, cy, ux, uy, r, t, gPx[j], gPy[j]))
                 count++;
         }
     }
     return  ( float )count / (N * M);
}

我先用伪随机数生成一些测试数据，然后用两个循环共执行1亿次检测。这比较合乎游戏应用的情况，每次通常是检测M个角色是否在1个扇形之内，做N个迭代。

使用 VS2010，缺省设置加上/fp:fast，Win32/Release的结果:

1
2
3
4
5
6
7
NoOperation hit=0% time =0.376s
IsPointInCircularSector hit=30.531% time =3.964s
 
NoOperation hit=0% time =0.364s
IsPointInCircularSector1 hit=30.531% time =0.643s
IsPointInCircularSector2 hit=30.531% time =0.614s
IsPointInCircularSector3 hit=30.531% time =0.571s

NoOperation是指测试的函数只简单传回false，用来检测函循环、读取测试数据和函数调用的开销。之后测的时间会减去这个开销。

可以看到，在简单优化后，预计算  和 ，并延后sqrt()，性能大幅提升4倍以上。之后的优化(第2和第3个版本)，只能再优化少许。为了减少sqrt()，却增加了乘数和比较。

有趣的是，若打开/arch:SSE2

1
2
3
4
5
6
7
NoOperation hit=0% time =0.453s
IsPointInCircularSector hit=30.531% time =2.645s
 
NoOperation hit=0% time =0.401s
IsPointInCircularSector1 hit=30.531% time =0.29s
IsPointInCircularSector2 hit=30.531% time =0.32s
IsPointInCircularSector3 hit=30.531% time =0.455s

所有性能都提升了，但优化版本反而一直倒退。主要原因应该是SSE含sqrtss指令，能快速完成开平方运算，第2个和第3个版本所做的优化反而增加了指令及分支，而第3个版本更需要把SSE寄存器的值储存至普通寄在器做整数运算，造成更大损失。

优化版本4和5: SOA SIMD

由于编译器自动化的SSE标量可以提高性，进一步的，可采用SOA(struct of array)布局进行以4组参数为单位的SIMD版本。使用了SSE/SSE2 intrinsic的实现，一个版本使用_mm_sqrt_ps()，一个版本去除了_mm_sqrt_ps():

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
// SSE2, SOA(struct of array) layout
// SSE2, SOA(struct of array) layout
__m128 IsPointInCircularSector4(
     __m128 cx, __m128 cy, __m128 ux, const  __m128& uy, const  __m128& squaredR, const  __m128& cosTheta,
     const  __m128& px, const  __m128& py)
{
     // D = P - C
     __m128 dx = _mm_sub_ps(px, cx);
     __m128 dy = _mm_sub_ps(py, cy);
 
     // |D|^2 = (dx^2 + dy^2)
     __m128 squaredLength = _mm_add_ps(_mm_mul_ps(dx, dx), _mm_mul_ps(dy, dy));
 
     // |D|^2 < r^2
     __m128 lengthResult = _mm_cmpgt_ps(squaredR, squaredLength);
 
     // |D|
     __m128 length = _mm_sqrt_ps(squaredLength);
 
     // D dot U
     __m128 DdotU = _mm_add_ps(_mm_mul_ps(dx, ux), _mm_mul_ps(dy, uy));
 
     // D dot U > |D| cos(theta)
     __m128 angularResult = _mm_cmpgt_ps(DdotU, _mm_mul_ps(length, cosTheta));
 
     __m128 result = _mm_and_ps(lengthResult, angularResult);
 
     return  result;
}
 
// SSE2, SOA(struct of array) layout without sqrt()
__m128 IsPointInCircularSector5(
     __m128 cx, __m128 cy, __m128 ux, const  __m128& uy, const  __m128& squaredR, const  __m128& cosTheta,
     const  __m128& px, const  __m128& py)
{
     // D = P - C
     __m128 dx = _mm_sub_ps(px, cx);
     __m128 dy = _mm_sub_ps(py, cy);
 
     // |D|^2 = (dx^2 + dy^2)
     __m128 squaredLength = _mm_add_ps(_mm_mul_ps(dx, dx), _mm_mul_ps(dy, dy));
 
     // |D|^2 < r^2
     __m128 lengthResult = _mm_cmpgt_ps(squaredR, squaredLength);
 
     // D dot U
     __m128 DdotU = _mm_add_ps(_mm_mul_ps(dx, ux), _mm_mul_ps(dy, uy));
 
     // D dot U > |D| cos(theta)
     // <=>
     // (D dot U)^2 > |D|^2 (cos(theta))^2 if D dot U >= 0 and cos(theta) >= 0
     // (D dot U)^2 < |D|^2 (cos(theta))^2 if D dot U <  0 and cos(theta) <  0
     // true                               if D dot U >= 0 and cos(theta) <  0
     // false                              if D dot U <  0 and cos(theta) >= 0
     __m128 asign = _mm_and_ps(DdotU, cSignMask.v);
     __m128 bsign = _mm_and_ps(cosTheta, cSignMask.v);
     __m128 equalsign = _mm_castsi128_ps(_mm_cmpeq_epi32(_mm_castps_si128(asign), _mm_castps_si128(bsign)));
     
     __m128 lhs = _mm_or_ps(_mm_mul_ps(DdotU, DdotU), asign);
     __m128 rhs = _mm_or_ps(_mm_mul_ps(squaredLength, _mm_mul_ps(cosTheta, cosTheta)), asign);
     __m128 equalSignResult = _mm_cmpgt_ps(lhs, rhs);
 
     __m128 unequalSignResult = _mm_castsi128_ps(_mm_cmpeq_epi32(_mm_castps_si128(asign), _mm_setzero_si128()));
 
     __m128 result = _mm_and_ps(lengthResult, _mm_or_ps(
         _mm_and_ps(equalsign, equalSignResult),
         _mm_andnot_ps(equalsign, unequalSignResult)));
 
     return  result;
}

因为要以SIMD并行执行，不能分支。所有分支都要计算，然后再混合(blend)在一起。测试结果:

1
2
IsPointInCircularSector4 hit=30.531% time =0.121s
IsPointInCircularSector5 hit=30.531% time =0.177s

SOA的好处是移植简单，基本上只要把float改成__m128，然后用相应的intrinsic代替浮点四则运算便可。而且能尽用每指令4个数据的吞吐量(throughput)。如果用__m128表示二维矢量(或其他问题中的三维矢量)，许多时候会浪费了一些吞吐量(例如在点积运算上)。

第4版本大约比之前最快的版本快近两倍，比没有/arch:SSE2最慢的版本快32倍。第5版本虽去除了开平方，在这个架构下还是得不偿失。

如果计上NoOperationSIMD的时间，SOA、对齐的内存访问方式，以及减少迭代数目，性能提高的幅度就更大。在实际应用中，通常比较难直套用这种内存布局，可能需要做转置(transpose)。对于开发高性能的子系统，在数据结构上可考虑直接储存为SOA。

其他可行方法

除了本文所述的方法，还可以用其他方法:

  • 把  转换至扇形的局部坐标  (如  为局部坐标的x轴)，然后再比较角度。这个是原来同事的做法。
  • 利用二维的"叉积"(实际上叉积只定义于三维)去判断  是否在两边之内。但需要把u旋转\pm \theta来计算两边。另外，边与  的夹角可能大于  ，要小心处理。
  • 查表。对于现在的机器来说，许多时候直接计算会比查表快，对于SOA来说，查表不能并行，更容易做成瓶颈。此外，查表要结合逻辑所需的精度，通用性会较弱。然言，对于某些应用场景，例如更复杂的运算、非常受限的处理单元(CPU/GPU/DSP等)，查表还是十分有用的思路。

或许你还会想到其他方法。也可尝试用不同的扇形表示方式。

结语

这个看来非常简单的问题，也可以作不同尝试。在日常的工作中，可以改善的地方就更是数之不尽了。


#include <cassert>
  #define _USE_MATH_DEFINES
  #include <cmath>
  #include <iostream>
  #include <ctime>
  #include <emmintrin.h>
   
  // Naive
  bool IsPointInCircularSector(
  float cx, float cy, float ux, float uy, float r, float theta,
  float px, float py)
  {
  assert(theta > 0 && theta < M_PI);
  assert(r > 0.0f);
   
  // D = P - C
  float dx = px - cx;
  float dy = py - cy;
   
  // |D| = (dx^2 + dy^2)^0.5
  float length = sqrt(dx * dx + dy * dy);
   
  // |D| > r
  if (length > r)
  return false;
   
  // Normalize D
  dx /= length;
  dy /= length;
   
  // acos(D dot U) < theta
  return acos(dx * ux + dy * uy) < theta;
  }
   
  // Basic: use squareR and cosTheta as parameters, defer sqrt(), eliminate division
  bool IsPointInCircularSector1(
  float cx, float cy, float ux, float uy, float squaredR, float cosTheta,
  float px, float py)
  {
  assert(cosTheta > -1 && cosTheta < 1);
  assert(squaredR > 0.0f);
   
  // D = P - C
  float dx = px - cx;
  float dy = py - cy;
   
  // |D|^2 = (dx^2 + dy^2)
  float squaredLength = dx * dx + dy * dy;
   
  // |D|^2 > r^2
  if (squaredLength > squaredR)
  return false;
   
  // |D|
  float length = sqrt(squaredLength);
   
  // D dot U > |D| cos(theta)
  return dx * ux + dy * uy > length * cosTheta;
  }
   
  // Eliminate sqrt()
  bool IsPointInCircularSector2(
  float cx, float cy, float ux, float uy, float squaredR, float cosTheta,
  float px, float py)
  {
  assert(cosTheta > -1 && cosTheta < 1);
  assert(squaredR > 0.0f);
   
  // D = P - C
  float dx = px - cx;
  float dy = py - cy;
   
  // |D|^2 = (dx^2 + dy^2)
  float squaredLength = dx * dx + dy * dy;
   
  // |D|^2 > r^2
  if (squaredLength > squaredR)
  return false;
   
  // D dot U
  float DdotU = dx * ux + dy * uy;
   
  // D dot U > |D| cos(theta)
  // <=>
  // (D dot U)^2 > |D|^2 (cos(theta))^2 if D dot U >= 0 and cos(theta) >= 0
  // (D dot U)^2 < |D|^2 (cos(theta))^2 if D dot U < 0 and cos(theta) < 0
  // true if D dot U >= 0 and cos(theta) < 0
  // false if D dot U < 0 and cos(theta) >= 0
  if (DdotU >= 0 && cosTheta >= 0)
  return DdotU * DdotU > squaredLength * cosTheta * cosTheta;
  else if (DdotU < 0 && cosTheta < 0)
  return DdotU * DdotU < squaredLength * cosTheta * cosTheta;
  else
  return DdotU >= 0;
  }
   
  // Bit trick
  bool IsPointInCircularSector3(
  float cx, float cy, float ux, float uy, float squaredR, float cosTheta,
  float px, float py)
  {
  assert(cosTheta > -1 && cosTheta < 1);
  assert(squaredR > 0.0f);
   
  // D = P - C
  float dx = px - cx;
  float dy = py - cy;
   
  // |D|^2 = (dx^2 + dy^2)
  float squaredLength = dx * dx + dy * dy;
   
  // |D|^2 > r^2
  if (squaredLength > squaredR)
  return false;
   
  // D dot U
  float DdotU = dx * ux + dy * uy;
   
  // D dot U > |D| cos(theta)
  // <=>
  // (D dot U)^2 > |D|^2 (cos(theta))^2 if D dot U >= 0 and cos(theta) >= 0
  // (D dot U)^2 < |D|^2 (cos(theta))^2 if D dot U < 0 and cos(theta) < 0
  // true if D dot U >= 0 and cos(theta) < 0
  // false if D dot U < 0 and cos(theta) >= 0
  const unsigned cSignMask = 0x80000000;
  union {
  float f;
  unsigned u;
  }a, b, lhs, rhs;
  a.f = DdotU;
  b.f = cosTheta;
  unsigned asign = a.u & cSignMask;
  unsigned bsign = b.u & cSignMask;
  if (asign == bsign) {
  lhs.f = DdotU * DdotU;
  rhs.f = squaredLength * cosTheta * cosTheta;
  lhs.u |= asign;
  rhs.u |= asign;
  return lhs.f > rhs.f;
  }
  else
  return asign == 0;
  }
   
  const union VectorI32 {
  unsigned u[4];
  __m128 v;
  }cSignMask = { 0x80000000, 0x80000000, 0x80000000, 0x80000000 };
   
  // SSE2, SOA(struct of array) layout
  __m128 IsPointInCircularSector4(
  __m128 cx, __m128 cy, __m128 ux, const __m128& uy, const __m128& squaredR, const __m128& cosTheta,
  const __m128& px, const __m128& py)
  {
  // D = P - C
  __m128 dx = _mm_sub_ps(px, cx);
  __m128 dy = _mm_sub_ps(py, cy);
   
  // |D|^2 = (dx^2 + dy^2)
  __m128 squaredLength = _mm_add_ps(_mm_mul_ps(dx, dx), _mm_mul_ps(dy, dy));
   
  // |D|^2 < r^2
  __m128 lengthResult = _mm_cmpgt_ps(squaredR, squaredLength);
   
  // |D|
  __m128 length = _mm_sqrt_ps(squaredLength);
   
  // D dot U
  __m128 DdotU = _mm_add_ps(_mm_mul_ps(dx, ux), _mm_mul_ps(dy, uy));
   
  // D dot U > |D| cos(theta)
  __m128 angularResult = _mm_cmpgt_ps(DdotU, _mm_mul_ps(length, cosTheta));
   
  __m128 result = _mm_and_ps(lengthResult, angularResult);
   
  return result;
  }
   
  // SSE2, SOA(struct of array) layout without sqrt()
  __m128 IsPointInCircularSector5(
  __m128 cx, __m128 cy, __m128 ux, const __m128& uy, const __m128& squaredR, const __m128& cosTheta,
  const __m128& px, const __m128& py)
  {
  // D = P - C
  __m128 dx = _mm_sub_ps(px, cx);
  __m128 dy = _mm_sub_ps(py, cy);
   
  // |D|^2 = (dx^2 + dy^2)
  __m128 squaredLength = _mm_add_ps(_mm_mul_ps(dx, dx), _mm_mul_ps(dy, dy));
   
  // |D|^2 < r^2
  __m128 lengthResult = _mm_cmpgt_ps(squaredR, squaredLength);
   
  // D dot U
  __m128 DdotU = _mm_add_ps(_mm_mul_ps(dx, ux), _mm_mul_ps(dy, uy));
   
  // D dot U > |D| cos(theta)
  // <=>
  // (D dot U)^2 > |D|^2 (cos(theta))^2 if D dot U >= 0 and cos(theta) >= 0
  // (D dot U)^2 < |D|^2 (cos(theta))^2 if D dot U < 0 and cos(theta) < 0
  // true if D dot U >= 0 and cos(theta) < 0
  // false if D dot U < 0 and cos(theta) >= 0
  __m128 asign = _mm_and_ps(DdotU, cSignMask.v);
  __m128 bsign = _mm_and_ps(cosTheta, cSignMask.v);
  __m128 equalsign = _mm_castsi128_ps(_mm_cmpeq_epi32(_mm_castps_si128(asign), _mm_castps_si128(bsign)));
   
  __m128 lhs = _mm_or_ps(_mm_mul_ps(DdotU, DdotU), asign);
  __m128 rhs = _mm_or_ps(_mm_mul_ps(squaredLength, _mm_mul_ps(cosTheta, cosTheta)), asign);
  __m128 equalSignResult = _mm_cmpgt_ps(lhs, rhs);
   
  __m128 unequalSignResult = _mm_castsi128_ps(_mm_cmpeq_epi32(_mm_castps_si128(asign), _mm_setzero_si128()));
   
  __m128 result = _mm_and_ps(lengthResult, _mm_or_ps(
  _mm_and_ps(equalsign, equalSignResult),
  _mm_andnot_ps(equalsign, unequalSignResult)));
   
  return result;
  }
   
  ///
  // Tests
   
  bool NoOperation(
  float cx, float cy, float ux, float uy, float squaredR, float cosTheta,
  float px, float py)
  {
  return false;
  }
   
  __m128 NoOperationSIMD(
  __m128 cx, __m128 cy, __m128 ux, const __m128& uy, const __m128& squaredR, const __m128& cosTheta,
  const __m128& px, const __m128& py)
  {
  return _mm_setzero_ps();
  }
   
  bool AreSame(
  float cx, float cy, float ux, float uy, float squaredR, float cosTheta,
  float px, float py)
  {
  bool result1 = IsPointInCircularSector1(
  cx, cy, ux, uy, squaredR, cosTheta, px, py);
   
  bool result2 = IsPointInCircularSector2(
  cx, cy, ux, uy, squaredR, cosTheta, px, py);
   
  bool result3 = IsPointInCircularSector3(
  cx, cy, ux, uy, squaredR, cosTheta, px, py);
   
  assert(result1 == result2);
  assert(result1 == result3);
  return result1;
  }
   
  void IsExpected(
  bool expected,
  float cx, float cy, float ux, float uy, float squaredR, float cosTheta,
  float px, float py)
  {
  assert(expected == AreSame(cx, cy, ux, uy, squaredR, cosTheta, px, py));
  }
   
  // Test data
   
  float Uniform(float minimum, float maximum) {
  return rand() * (maximum - minimum) / RAND_MAX + minimum;
  }
   
  const unsigned N = 1000;
  #if _DEBUG
  const unsigned M = 1000;
  #else
  const unsigned M = 100000;
  #endif
   
  __declspec(align(128)) float gCx[N];
  __declspec(align(128)) float gCy[N];
  __declspec(align(128)) float gUx[N];
  __declspec(align(128)) float gUy[N];
  __declspec(align(128)) float gSquaredR[N];
  __declspec(align(128)) float gCosTheta[N];
  __declspec(align(128)) float gPx[M];
  __declspec(align(128)) float gPy[M];
  __declspec(align(128)) float gR[N]; // naive
  __declspec(align(128)) float gTheta[N]; // naive
   
  // Run test
   
  template <bool naiveParam, typename TestFunc>
  float Test(TestFunc f, float rmax = 2.0f) {
  unsigned count = 0;
  for (int i = 0; i < N; i++) {
  float cx = gCx[i];
  float cy = gCy[i];
  float ux = gUx[i];
  float uy = gUy[i];
  float r = naiveParam ? gR[i] : gSquaredR[i];
  float t = naiveParam ? gTheta[i] : gCosTheta[i];
   
  for (int j = 0; j < M; j++) {
  if (f(cx, cy, ux, uy, r, t, gPx[j], gPy[j]))
  count++;
  }
  }
  return (float)count / (N * M);
  }
   
  template <typename TestFunc>
  float TestSIMD(TestFunc f, float rmax = 2.0f) {
  static const unsigned cCountTable[] = {
  0, // 0000
  1, // 0001
  1, // 0010
  2, // 0011
  1, // 0100
  2, // 0101
  2, // 0110
  3, // 0111
  1, // 1000
  2, // 1001
  2, // 1010
  3, // 1011
  2, // 1100
  3, // 1101
  3, // 1110
  4, // 1111
  };
   
  unsigned count = 0;
  for (int i = 0; i < N; i++) {
  __m128 cx = _mm_set1_ps(gCx[i]);
  __m128 cy = _mm_set1_ps(gCy[i]);
  __m128 ux = _mm_set1_ps(gUx[i]);
  __m128 uy = _mm_set1_ps(gUy[i]);
  __m128 squaredR = _mm_set1_ps(gSquaredR[i]);
  __m128 cosTheta = _mm_set1_ps(gCosTheta[i]);
   
  for (int j = 0; j < M; j += 4) {
  __m128 px = _mm_load_ps(&gPx[j]);
  __m128 py = _mm_load_ps(&gPy[j]);
  int mask = _mm_movemask_ps(f(cx, cy, ux, uy, squaredR, cosTheta, px, py));
  count += cCountTable[mask];
  }
  }
  return (float)count / (N * M);
  }
   
  template <bool naiveParam, typename TestFunc>
  float Performance(const char* name, TestFunc f, float overhead = 0.0f) {
  clock_t start = clock();
  float hit = Test<naiveParam>(f);
  clock_t end = clock();
  float time = (float)(end - start) / CLOCKS_PER_SEC - overhead;
  printf("%s hit=%g%% time=%gs\n", name, hit * 100.0f, time);
  return time;
  }
   
  template <typename TestFunc>
  float PerformanceSIMD(const char* name, TestFunc f, float overhead = 0.0f) {
  clock_t start = clock();
  float hit = TestSIMD(f);
  clock_t end = clock();
  float time = (float)(end - start) / CLOCKS_PER_SEC - overhead;
  printf("%s hit=%g%% time=%gs\n", name, hit * 100.0f, time);
  return time;
  }
   
  int main(int argc, char* argv[]) {
  srand(0);
  for (int i = 0; i < N; i++) {
  gCx[i] = Uniform(-1.0f, 1.0f);
  gCy[i] = Uniform(-1.0f, 1.0f);
   
  float ux = Uniform(-1.0f, 1.0f);
  float uy = Uniform(-1.0f, 1.0f);
  float ru = 1.0f / (sqrt(ux * ux + uy * uy));
  gUx[i] = ux * ru;
  gUy[i] = uy * ru;
   
  gR[i] = Uniform(0.0f, 2.0f);
  gSquaredR[i] = gR[i] * gR[i];
  gTheta[i] = Uniform(0.0f, float(M_PI));
  gCosTheta[i] = cos(gTheta[i]);
  }
   
  for (int j = 0; j < N; j++) {
  gPx[j] = Uniform(-1.0f, 1.0f);
  gPy[j] = Uniform(-1.0f, 1.0f);
  }
   
  #ifdef _DEBUG
  IsExpected(true, 0.0f, 0.0f, 1.0f, 0.0f, 4.0f, 0.5f, 1.0f, 0.0f);
  IsExpected(false, 0.0f, 0.0f, 1.0f, 0.0f, 4.0f, 0.5f, 4.0f, 1.0f);
  IsExpected(true, 0.0f, 0.0f, 1.0f, 0.0f, 1.0f, 0.0f, 0.5f, 0.0f);
  IsExpected(false, 0.0f, 0.0f, 1.0f, 0.0f, 1.0f, 0.0f, -0.5f, 0.0f);
   
  Test<false>(AreSame);
  #else
  #define PERF(f, naiveParam, overhead) Performance<naiveParam>(#f, f, overhead)
  #define PERFSIMD(f, overhead) PerformanceSIMD(#f, f, overhead)
  float overhead = PERF(NoOperation, true, 0.0f);
  PERF(IsPointInCircularSector, true, overhead);
  printf("\n");
   
  overhead = PERF(NoOperation, false, 0.0f);
  PERF(IsPointInCircularSector1, false, overhead);
  PERF(IsPointInCircularSector2, false, overhead);
  PERF(IsPointInCircularSector3, false, overhead);
  printf("\n");
   
  overhead = PERFSIMD(NoOperationSIMD, 0.0f);
  PERFSIMD(IsPointInCircularSector4, overhead);
  PERFSIMD(IsPointInCircularSector5, overhead);
  #undef PERF
  #undef PERFSIMD
  #endif
  return 0;
  }

警醒与鞭策
关注
  • 0
    点赞
  • 踩
  • 1
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
专栏目录
面试题：检测点是否在扇形之内
weixin_33695450的博客
04-19 839
前几天，同事在报告中提及检测角色是否在扇形攻击范围的方法。我觉得该方法的性能不是太好，提出另一个颇为直接的方法。此问题在游戏中十分常见，只涉及简单的数学，却又可以看出实现者是否细心，所以我觉得可当作一道简单的面试题。问题在微博发表后得到不少回应，故撰文提供一些解答。 问题定义: 在二维中，检测点\mathbf{p}是否在扇形(circular sector)内，设扇形的顶点为\math...
判断一个点是否落在扇形区域内
qq_35451217的博客
12-25 4006
工作需要判定一个坐标（x,y）是否落在扇形区域内； 使用到数学函数库 #include<math.h> atan2(y,x)输出为弧度 返回角度： double cacul_alpha(double x,double y) { return atan2(y,x)*180/PI; } 判断一个点是否落在扇形区域内 char judge_coor_pos_fan...
安卓判断某个点是否在扇形或者区域内
07-19
安卓判断某个点是否在扇形或者区域内，系统自带判断，效率高
判断点在扇形内
残獄の世界に救いを
05-13 2084
原理 假设已知扇形半径 弧度 朝向 圆心点 想判断任意点P是否在扇形内,需要两个条件: 1.点P与圆心O的距离不大于 朝向与扇形的交点A到圆心O的距离 2.角AOP的度数不大于扇形度数/2 代码 --- @brief 判断某个点是否在扇形内 --- @public --- @param _degree number 扇形的弧度 --- @param _radius number 扇形的半径 --- @param _startPos UnityEngine.Vector3 角色所在的位置 扇形的圆心 ---
判断点是否在扇形范围内
zgjllf1011的博客
03-09 1万+
最近做一个潜入游戏其中有一种敌人面前会有一个扇形的检测范围，并且此扇形范围会来回旋转移动，所以就遇到了一个如何判断主角是否在扇形范围内的问题。在网上找了一些方法，但是感觉都太复杂，所以这里自己实现了一下。思路如下：我的扇形的生成方式是从中心方向像两边展开若干角度。所以我们只要算出扇形的中心向量与角色向量的夹脚是否大于扇形的展开角度即可。步骤：扇形的来回摆动是有程序控制，所以这里可以很容易得到扇形当...
Unity实现图形相交检测
08-19
在 Unity 游戏引擎中，图形相交检测广泛应用于游戏开发中，如伤害判定、碰撞检测等。下面将详细介绍 Unity 实现图形相交检测的相关知识点。 圆形与圆形的相交检测 圆形与圆形的相交检测是最简单的一种图形相交检测...
2015_2016七年级数学上册4.4扇形统计图检测无答案新版青岛版
08-06
14. **扇形统计图的计算**：（1）利用扇形统计图的百分比之和为100%，可求得抽查的学生总数；（2）参加“音乐”活动项目的百分比乘以360度，可以得到对应的圆心角度数；（3）利用“美术”活动项目在扇形图中的百分比...
Unity3D实现攻击范围检测
08-19
在上面的代码中，我们首先计算玩家与敌人的距离，然后计算玩家与敌人的方向向量，最后通过点积的计算来判断敌人是否在攻击范围内。 半圆形攻击范围检测 半圆形攻击范围检测是 Unity3D 中另一种攻击范围检测方法。...
基于3D LiDAR的矿山无人驾驶车行驶边界检测
06-10
准确的路沿检测与路沿图构建同时也是高精地图构建中的首要步骤，通过机器学的方法可以在现有路沿点的基础上初步构建道路结构拓扑。与传统的卡尔曼框架下的路沿跟踪所不同的是，本文提出了基于占用栅格思想的路沿跟踪...
江苏省常州市武进区九年级数学上册2.7弧长及扇形的面积课堂学习检测题二新版苏科版
08-06
这些题目涉及的是初中数学中的弧长和扇形面积的相关知识点，主要涵盖圆锥侧面积的计算、角度变化对形状的影响、直角三角形和圆的综合应用、扇形面积的求解以及图形旋转与路径长度的计算。下面将对这些知识点进行详细...
判断点是否在扇形区域内
12-29
判断点是否在扇形区域内
android判断点击位置是否在扇形区域内
08-27
主要为大家详细介绍了android判断点击位置是否在扇形区域内，具有一定的参考价值，感兴趣的小伙伴们可以参考一下
android 扇形点击范围,android判断点击位置是否在扇形区域内
weixin_35554306的博客
05-26 245
在做仿支付宝记账本界面效果时遇到了一个问题，在环形图中点击每一个环是会显示出不同的内容，因此，必须判断用户到底点击了哪个圆环，网上查阅资料说可以根据颜色来判断，但是心里总是觉得根据颜色不是很好的解决方案(虽然可以解决这个问题)，因此记录下自己解决次问题的方法，也有助于自己后续复习。绘制扇形要判断点击的位置是否是在扇形区域内，首先要绘制扇形，绘制扇形的方法可以使用canvas.drawArc()方法...
判断点是否在扇形内解法探究 C++实现
qq_40515692的博客
09-05 1851
#include <iostream> using namespace std; const int n = 50,m = 50; bool table[n][m]; void func(int x,int y,int r){ for(int ex=0;ex<n;ex++) for(int ey=0;ey<m;ey++){ int dx = x-ex; int dy = y-ey; if(dx*dx + dy*dy > r*r){
Unity检测目标点是否在扇形攻击范围内
热门推荐
改变世界的力量
07-09 1万+
在游戏设计中,有的技能的攻击范围是一个扇形,判断场景中的怪物是否在扇形范围内,我们可以理解具体思路: 首先,判断该目标点与技能释放者的距离是否小于扇形的半径; 其次目标点与技能释放者面对的方向向量夹角是否小于扇形角度的二分之一; 如上图所示,1为技能释放者,2为目标点 代码如下: using UnityEngine; using System.Collections
three.js点乘.dot()判断是否在扇形内
u014291990的博客
01-01 875
在扇形范围内，从人所在位置，向外沿着扇形半径方向绘制向量c，你会发现向量c与向量a最大夹角就是扇形角度一半。上节课给大家总结过，人指向物体的向量b和表示人正前方的向量a加载在0~180度之间。换句话说，向量b与向量a夹角越小，对应的余弦值越大。比较距离L与扇形半径的关系，排除物体不在扇形范围内的情况。角度θ在0 ~ 180度范围内，随着θ度数的增加，角度θ。判断一个点(物体)是否在人前方扇形范围内(半径。判断一个点(物体)是否在人前方的扇形范围内。点乘计算向量a和b的夹角余弦值cos。
java 点在扇形内_面试题：检测点是否在扇形之内
weixin_39892311的博客
02-27 696
前几天，同事在报告中提及检测角色是否在扇形攻击范围的方法。我觉得该方法的性能不是太好，提出另一个颇为直接的方法。此问题在游戏中十分常见，只涉及简单的数学，却又可以看出实现者是否细心，所以我觉得可当作一道简单的面试题。问题在微博发表后得到不少回应，故撰文提供一些解答。问题定义:在二维中，检测点\mathbf{p}是否在扇形(circular sector)内，设扇形的顶点为\mathbf{c}，半径...
图形问题
qq_33666011的博客
03-31 557
转转：不记得从哪里转载的了，如果博主发现，请私信告知 已知两圆圆心坐标和半径，求相交部分面积 //如果两个圆的圆心距离+小半径，小于 大半径，就说明两个圆内交或包含，那么相交面积就是小圆面积 //如果两个圆的圆心距离>大半径+小半径，说明两者外交或相离，那么相交面积就是0 //前两者都不是说明，两者相交。相交的面积等于，（r1的扇形-三角形）+（r2的扇形-三角形） //求三角形...
unity3d 射线检测 扇形
最新发布
01-25
而射线检测扇形则是在一定范围内进行射线检测，常用于判断某个角度内是否有其他物体存在。 首先，我们需要确定扇形的起始点和扇形的角度范围。然后，通过编写代码来创建射线，并设置射线的起点和方向，以及射线的...
评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值