『代码』手写IOU计算
计算2d iou
来自OpenPCDet的eval.py。关键就在于：重叠部分的宽度 = 两个图像右边缘的较小值 - 两个图像左边缘的较大值，重叠部分的高度 = 两个图像下边缘的较小值 - 两个图像上边缘的较大值U。右小减左大，下小减上大。然后IoU等于交除以并减交
点击查看代码
@numba.jit(nopython=True)
def image_box_overlap(boxes, query_boxes, criterion=-1):
"""
计算图像box的iou
Args:
boxes:一个part中的全部gt，以第一个part为例(642,4)
query_boxes：一个part中的全部dt，以第一个part为例(233,4)
"""
N = boxes.shape[0] # gt_box的总数
K = query_boxes.shape[0] # det_box的总数
# 初始化overlap矩阵
overlaps = np.zeros((N, K), dtype=boxes.dtype)
# 两层for循环逐个box计算iou，因为有jit加速，所以成for循环的形式
for k in range(K):
# 计算第k个dt box的面积（box是左上角和右下角的形式[x1,y1,x2,y2]）
qbox_area = ((query_boxes[k, 2] - query_boxes[k, 0]) *
(query_boxes[k, 3] - query_boxes[k, 1]))
for n in range(N): # 遍历gt boxes
# 重叠部分的宽度 = 两个图像右边缘的较小值 - 两个图像左边缘的较大值
iw = (min(boxes[n, 2], query_boxes[k, 2]) -
max(boxes[n, 0], query_boxes[k, 0]))
if iw > 0: # 如果宽度方向有重叠，再计算高度
# 重叠部分的高度 = 两个图像下边缘的较小值 - 两个图像上边缘的较大值
ih = (min(boxes[n, 3], query_boxes[k, 3]) -
max(boxes[n, 1], query_boxes[k, 1]))
if ih > 0:
if criterion == -1: # 默认执行criterion = -1
# 求两个box的并
ua = (
(boxes[n, 2] - boxes[n, 0]) *
(boxes[n, 3] - boxes[n, 1]) + qbox_area - iw * ih)
elif criterion == 0:
ua = ((boxes[n, 2] - boxes[n, 0]) *
(boxes[n, 3] - boxes[n, 1]))
elif criterion == 1:
ua = qbox_area
else:
ua = 1.0
overlaps[n, k] = iw * ih / ua
return overlaps
计算3d iou
首先，结果可以拆分成2d bev下带rotation的IoU与height交叉高度相乘。于是，问题的关键其实在于2d rotated iou这件事。计算2d rotated iou：
1. 根据box信息计算corners（基于box center，先得到local的角点坐标，并将corners编码为顺时针顺序，然后旋转加平移得到绝对坐标）
2. 求两box范围的所有交点，于是得到polygon。交点虽然不规则，但总不会超过八个点，于是可以预先分配buffer。分为a. 在对方box内的点 b. 边与边的交叉点
3. 根据polygon中心点对角点进行排序，使其满足逆时针方向
4. 计算多边形面积。采用累计计算三角形面积（利用叉乘）的方式（由于夹角逐渐增大，与第一个点构成的三角形面积一定是在同侧的符号相同）
计算box的3d iou
从d3_box_overlap()开始，其包括： rotate_iou_gpu_eval() 计算bev iou + d3_box_overlap_kernel() 结合bev iou计算3d iou
点击查看代码
def d3_box_overlap(boxes, qboxes, criterion=-1):
"""
Args:
boxes:一个part中的全部gt，以第一个part为例(642,7) [x,y,z,dx,dy,dz,alpha]
query_boxes：一个part中的全部dt，以第一个part为例(233,7)
rinc:在鸟瞰图视角下的iou（642，233）
Returns:
返回3D iou-->rinc
"""
rinc = rotate_iou_gpu_eval(boxes[:, [0, 2, 3, 5, 6]],
qboxes[:, [0, 2, 3, 5, 6]], 2)
d3_box_overlap_kernel(boxes, qboxes, rinc, criterion)
return rinc
点击查看代码
@numba.jit(nopython=True, parallel=True)
def d3_box_overlap_kernel(boxes, qboxes, rinc, criterion=-1):
"""
计算box的3D iou
Args:
boxes:一个part中的全部gt，以第一个part为例(642,7) [x,y,z,dx,dy,dz,alpha]
query_boxes：一个part中的全部dt，以第一个part为例(233,7)
rinc:在鸟瞰图视角下的iou（642，233）
Returns:
返回3D iou-->rinc
"""
# ONLY support overlap in CAMERA, not lidar.
# 在相机坐标系下进行计算，z轴朝前，y轴朝下，x轴朝右，因此，俯视图取x和z
N, K = boxes.shape[0], qboxes.shape[0]
# 遍历gt
for i in range(N):
# 遍历dt
for j in range(K):
# 如果鸟瞰视角存在重叠
if rinc[i, j] > 0:
# iw = (min(boxes[i, 1] + boxes[i, 4], qboxes[j, 1] +
# qboxes[j, 4]) - max(boxes[i, 1], qboxes[j, 1]))
# 这里的1是y轴，在相机坐标系下就是高度方向，重叠高度=上边缘的最小值-下边缘的最大值
iw = (min(boxes[i, 1], qboxes[j, 1]) - max(
boxes[i, 1] - boxes[i, 4], qboxes[j, 1] - qboxes[j, 4]))
# 如果重叠高度 > 0
if iw > 0:
# 1.求两个box的体积
area1 = boxes[i, 3] * boxes[i, 4] * boxes[i, 5]
area2 = qboxes[j, 3] * qboxes[j, 4] * qboxes[j, 5]
# 2.求交集体积
inc = iw * rinc[i, j]
# 3.根据criterion，计算并集体积
if criterion == -1:
ua = (area1 + area2 - inc)
elif criterion == 0:
ua = area1
elif criterion == 1:
ua = area2
else:
ua = inc
# 4.计算交并比
rinc[i, j] = inc / ua
else:
rinc[i, j] = 0.0
计算2d bev iou
rotate_iou_gpu_eval() → rotate_iou_kernel_eval() → devRotateIoUEval()
@cuda.jit('(float32[:], float32[:], int32)', device=True, inline=True)
def devRotateIoUEval(rbox1, rbox2, criterion=-1):
# 1.计算两个box的面积
area1 = rbox1[2] * rbox1[3]
area2 = rbox2[2] * rbox2[3]
# 2.求两个box的交集
area_inter = inter(rbox1, rbox2)
# 3.求交并比
if criterion == -1:
return area_inter / (area1 + area2 - area_inter)
elif criterion == 0:
return area_inter / area1
elif criterion == 1:
return area_inter / area2
else:
return area_inter
开始计算交集
点击查看代码
@cuda.jit('(float32[:], float32[:])', device=True, inline=True)
def inter(rbbox1, rbbox2):
# 0.初始化
# 在cuda上开辟局部内存，存储角点信息
corners1 = cuda.local.array((8, ), dtype=numba.float32)
corners2 = cuda.local.array((8, ), dtype=numba.float32)
# 由于交集部分不规则，但是最多8个点，共16个坐标
intersection_corners = cuda.local.array((16, ), dtype=numba.float32)
# 1.求两个box的四个角点坐标，先平移，然后旋转
rbbox_to_corners(corners1, rbbox1)
rbbox_to_corners(corners2, rbbox2)
# 2.求box的所有交点
num_intersection = quadrilateral_intersection(corners1, corners2,
intersection_corners)
# 3.根据对交集多边形中心点对其角点进行排序，使其满足逆时针方向
sort_vertex_in_convex_polygon(intersection_corners, num_intersection)
# print(intersection_corners.reshape([-1, 2])[:num_intersection])
# 4.根据叉乘几何意义，求取交集多边形面积（由于夹角逐渐增大，与第一个点构成的三角形面积一定是在同侧的符号相同）
return area(intersection_corners, num_intersection)
1. 计算角点坐标
@cuda.jit('(float32[:], float32[:])', device=True, inline=True)
def rbbox_to_corners(corners, rbbox):
# generate clockwise corners and rotate it clockwise
# 将中心点加宽高的形式和角度的形式转换为角点形式(x1,y1,x2,y2,x3,y3,x4,y4)
# 1.初始化
angle = rbbox[4]
a_cos = math.cos(angle)
a_sin = math.sin(angle)
center_x = rbbox[0]
center_y = rbbox[1]
x_d = rbbox[2]
y_d = rbbox[3]
corners_x = cuda.local.array((4, ), dtype=numba.float32)
corners_y = cuda.local.array((4, ), dtype=numba.float32)
# 2.计算在坐标原点的box模板
corners_x[0] = -x_d / 2
corners_x[1] = -x_d / 2
corners_x[2] = x_d / 2
corners_x[3] = x_d / 2
corners_y[0] = -y_d / 2
corners_y[1] = y_d / 2
corners_y[2] = y_d / 2
corners_y[3] = -y_d / 2
# 3.然后旋转和平移，计算4个角点
for i in range(4):
corners[2 *
i] = a_cos * corners_x[i] + a_sin * corners_y[i] + center_x
corners[2 * i
+ 1] = -a_sin * corners_x[i] + a_cos * corners_y[i] + center_y
2. 计算交集多边形：在框内的点 + 边与边的交点
@cuda.jit('(float32[:], float32[:], float32[:])', device=True, inline=True)
def quadrilateral_intersection(pts1, pts2, int_pts):
"""
计算两个box交集的角点坐标
Args:
pts1:(8,)
pts2:(8,)
Returns:
int_pts:(16,)
"""
num_of_inter = 0
# 计算交集的角点
for i in range(4):
# 分别检查四个角点在对方box内
if point_in_quadrilateral(pts1[2 * i], pts1[2 * i + 1], pts2):
int_pts[num_of_inter * 2] = pts1[2 * i]
int_pts[num_of_inter * 2 + 1] = pts1[2 * i + 1]
num_of_inter += 1
if point_in_quadrilateral(pts2[2 * i], pts2[2 * i + 1], pts1):
int_pts[num_of_inter * 2] = pts2[2 * i]
int_pts[num_of_inter * 2 + 1] = pts2[2 * i + 1]
num_of_inter += 1
temp_pts = cuda.local.array((2, ), dtype=numba.float32)
for i in range(4):
for j in range(4):
# 求交点坐标
has_pts = line_segment_intersection(pts1, pts2, i, j, temp_pts)
if has_pts:
int_pts[num_of_inter * 2] = temp_pts[0]
int_pts[num_of_inter * 2 + 1] = temp_pts[1]
num_of_inter += 1
return num_of_inter
3. 排序多边形角点
@cuda.jit('(float32[:], int32)', device=True, inline=True)
def sort_vertex_in_convex_polygon(int_pts, num_of_inter):
if num_of_inter > 0:
center = cuda.local.array((2, ), dtype=numba.float32)
center[:] = 0.0
# 1.计算多边形的平均中心点
for i in range(num_of_inter):
center[0] += int_pts[2 * i]
center[1] += int_pts[2 * i + 1]
center[0] /= num_of_inter
center[1] /= num_of_inter
# 2.对多边形的点排序
v = cuda.local.array((2, ), dtype=numba.float32)
vs = cuda.local.array((16, ), dtype=numba.float32)
for i in range(num_of_inter):
v[0] = int_pts[2 * i] - center[0]
v[1] = int_pts[2 * i + 1] - center[1]
d = math.sqrt(v[0] * v[0] + v[1] * v[1])
v[0] = v[0] / d
v[1] = v[1] / d
if v[1] < 0:
v[0] = -2 - v[0]
vs[i] = v[0]
j = 0
temp = 0
for i in range(1, num_of_inter):
if vs[i - 1] > vs[i]:
temp = vs[i]
tx = int_pts[2 * i]
ty = int_pts[2 * i + 1]
j = i
while j > 0 and vs[j - 1] > temp:
vs[j] = vs[j - 1]
int_pts[j * 2] = int_pts[j * 2 - 2]
int_pts[j * 2 + 1] = int_pts[j * 2 - 1]
j -= 1
vs[j] = temp
int_pts[j * 2] = tx
int_pts[j * 2 + 1] = ty
4. 计算多边形面积
@cuda.jit('(float32[:], int32)', device=True, inline=True)
def area(int_pts, num_of_inter):
area_val = 0.0
# 逐多边形点遍历，求解面积
for i in range(num_of_inter - 2):
area_val += abs(
trangle_area(int_pts[:2], int_pts[2 * i + 2:2 * i + 4],
int_pts[2 * i + 4:2 * i + 6]))
return area_val
@cuda.jit('(float32[:], float32[:], float32[:])', device=True, inline=True)
def trangle_area(a, b, c):
"""
计算三角形面积（叉乘法）
/a
c/___b
"""
return ((a[0] - c[0]) * (b[1] - c[1]) - (a[1] - c[1]) *(b[0] - c[0])) / 2.0
浙公网安备 33010602011771号