2021-11-04

Ray-AABB交叉检测算法

——判断射线是否穿过立方体AABB

slab的碰撞检测算法(slab method)

---

1. 定义

AABB：Axially Aligned Bounding Box（轴对齐矩形边界框）

slab(平板)：在此指两个平行平面或直线之间的空间。3D空间中则是两坐标平面(xy面/xz面/yz面)之间的区域。由此，3D中的AABB盒子看作是由三组平行面形成的三个方向的Slab交集。

候选面：3D空间中，一组平行面正对射线的面叫候选面，另一面叫背面

---

2. 性质

• 如果一个点在AABB中，那么这个点必定同时在这3个slab中。

• 如果一条射线和AABB相交，那么这条射线和3个slab的相交部分必定有重合部分。

• 当射线与这三个候选面中的一个发生交叉之后，射线Ray的原点到这个面的距离要比到其他几个面的距离要长。

从图中我们可以得到：
(1) t2>t1
（其中t表示射线原点到交叉点的距离）
(2) max(t1,t2)<=tA<=min(t3,t4)

将2D推广到3D空间，有
max(t1,t2,t3)<=min(t4,t5,t6)
如果发生交叉，射线原点在射线方向到最近交叉面的距离为max(t1,t2,t3)

---

3. 步骤
   确定A射线与AABB是否交叉的三步骤：
   （1） 确定候选面
   将平面方程带入射线Ray的方程，求出这两个平面的t值，然后t值较小的那个自然先与射线交叉，那么就表示它是一个候选面。射线可以用参数方程表示为R(t) = P0 + t·d, (其中P0为射线起点，d为射线的方向向量)
   （2） 确定候选面方程
   平面由隐式定义方程X·n=D
   (其中X为平面上的点，n为平面法向量，D为原点到平面的距离)
   (可由平面方程的点法式推出)
   （3） 判断交叉点是否在AABB盒上
   max(t1,t2,t3)<=min(t4,t5,t6)

4. 碰撞检测算法公式
   (1)射线方程： R(t) = P0+t*d (1)
   (2)平面方程： X *n = D (2)
   P0 : 射线起点
   d : 射线的方向向量
   X ：平面上的点
   n : 法向量
   D ： 坐标原点到平面距离
   射线与平面相交时，将(1)（2)联立，有
   $\left(P_0+t\cdot d\right)\cdot n=D\Rightarrow t=\frac{\left(D-P_0\cdot n\right)}{d\cdot n}$
   由于AABB的slab平行于坐标轴，公式化简，设P0=（Px,Py,Pz)，d=(dx,dy,dz)时

• t和x-slab面的交点:
  法向量(1,0,0)
  $t=\frac{D-P_x}{dx}$
• t和y-slab面的交点:
  法向量(0,1,0)
  $t=\frac{D-P_y}{dy}$
• t和z-slab面的交点:
  法向量法向量(0,0,1)
  $t=\frac{D-P_z}{dz}$

5. 代码
   (以下代码来自论文differential volumetric rending:learning implicit 3D representations without 3D supervision.

def check_ray_intersection_with_unit_cube(ray0, ray_direction, padding=0.1,
                                          eps=1e-6):
    ''' Checks if rays ray0 + d * ray_direction intersect with unit cube with
    padding padding.
    
    ## create a unit cube ,detect weather the ray is intersect with the unit cube

    It returns the two intersection points as well as the sorted ray lengths d.

    Args:
        ray0 (tensor): Start positions of the rays
        ray_direction (tensor): Directions of rays
        padding (float): Padding which is applied to the unit cube
        eps (float): The epsilon value for numerical stability
    '''
    batch_size, n_pts, _ = ray0.shape
    device = ray0.device

    # calculate intersections with unit cube (< . , . >  is the dot product)
    # <n, x - p> = <n, ray0 + d * ray_direction - p_e> = 0
    # d = - <n, ray0 - p_e> / <n, ray_direction>

    # Get points on plane p_e
    p_distance = 0.5 + padding/2 # 一半边长加上padding的一半
    p_e = torch.ones(batch_size, n_pts, 6).to(device) * p_distance # Unit_cube的六个面
    p_e[:, :, 3:] *= -1.
	#p_e表示原点到unit cube的距离
	
    # Calculate the intersection points with given formula
    nominator = p_e - ray0.repeat(1, 1, 2)
    denominator = ray_direction.repeat(1, 1, 2)
    d_intersect = nominator / denominator
    p_intersect = ray0.unsqueeze(-2) + d_intersect.unsqueeze(-1) * \
        ray_direction.unsqueeze(-2)
	由于将六个面都融合进一个矩阵，故在表示坐标的维度前增加一个表示面的维度
    # Calculate mask where points intersect unit cube
    p_mask_inside_cube = (
        (p_intersect[:, :, :, 0] <= p_distance + eps) &
        (p_intersect[:, :, :, 1] <= p_distance + eps) &
        (p_intersect[:, :, :, 2] <= p_distance + eps) &
        (p_intersect[:, :, :, 0] >= -(p_distance + eps)) &
        (p_intersect[:, :, :, 1] >= -(p_distance + eps)) &
        (p_intersect[:, :, :, 2] >= -(p_distance + eps))
    ).cpu()
	由边界条件判断射线是否和面intersect
	
    # Correct rays are these which intersect exactly 2 times
    mask_inside_cube = p_mask_inside_cube.sum(-1) == 2
	#正确intersect的射线应该有穿入和穿出Unit cube
	
    # Get interval values for p's which are valid
    p_intervals = p_intersect[mask_inside_cube][p_mask_inside_cube[
        mask_inside_cube]].view(-1, 2, 3)
    p_intervals_batch = torch.zeros(batch_size, n_pts, 2, 3).to(device)
    p_intervals_batch[mask_inside_cube] = p_intervals

    # Calculate ray lengths for the interval points
    d_intervals_batch = torch.zeros(batch_size, n_pts, 2).to(device)
    norm_ray = torch.norm(ray_direction[mask_inside_cube], dim=-1)
    d_intervals_batch[mask_inside_cube] = torch.stack([
        torch.norm(p_intervals[:, 0] -
                   ray0[mask_inside_cube], dim=-1) / norm_ray,
        torch.norm(p_intervals[:, 1] -
                   ray0[mask_inside_cube], dim=-1) / norm_ray,
    ], dim=-1)

    # Sort the ray lengths
    d_intervals_batch, indices_sort = d_intervals_batch.sort()
    p_intervals_batch = p_intervals_batch[
        torch.arange(batch_size).view(-1, 1, 1),
        torch.arange(n_pts).view(1, -1, 1),
        indices_sort
    ]

    return p_intervals_batch, d_intervals_batch, mask_inside_cube

本文参考了 Ray-AABB,若不慎侵权，请通知删除。