等矩形全景图转透视图像的畸变问题排查

我有一张比例2:1、分辨率4000×2000的等矩形全景图，尝试根据相机朝向向量、目标图像尺寸和2D视场角（FOV）将其转换为透视图像。但发现当朝向向量同时包含x和y分量时，生成的透视图像会出现旋转错误与畸变，其余朝向的转换结果均正常。相关代码如下：

import cupy as cp
from scipy import ndimage

def calculate_R(orientation, up=cp.array([0, -1, 0])):

    orientation = orientation / cp.linalg.norm(orientation)

    if cp.allclose(orientation, up):
        return cp.array([[1, 0, 0], [0, 0, 1], [0, -1, 0]])

    if cp.allclose(orientation, -up):
        return cp.array([[1, 0, 0], [0, 0, -1], [0, 1, 0]])

    x = cp.cross(up, orientation)
    x = x / cp.linalg.norm(x)

    y = cp.cross(orientation, x)
    y = y / cp.linalg.norm(y)

    R = cp.array([x, y, orientation])

    return R

def rotate_xyz(R, x, y, z):

    xyz = cp.stack((x.ravel(), y.ravel(), z.ravel()))
    x, y, z = cp.matmul(R, xyz).reshape(3, *x.shape)
    
    return x, y, z

def transform_back(pano, orientation, image_shape, fov=[90, 90]):

    orientation = orientation / cp.linalg.norm(orientation)

    fov_x = cp.radians(fov[0])
    fov_y = cp.radians(fov[1])

    x_max = cp.tan(fov_x / 2)
    y_max = cp.tan(fov_y / 2)

    x = cp.linspace(-x_max, x_max, num=image_shape[1], dtype=cp.float32)
    y = cp.linspace(-y_max, y_max, num=image_shape[0], dtype=cp.float32)
    x, y = cp.meshgrid(x, y)
    z = -cp.ones_like(x)

    R = calculate_R(orientation)

    x, y, z = rotate_xyz(R, x, y, z)
    xyz = cp.stack([x, y, z], axis=-1).reshape(-1, 3)

    lon = cp.arctan2(xyz[:, 0], xyz[:, 2])
    lat = cp.arcsin(xyz[:, 1] / cp.linalg.norm(xyz, axis=1))

    u = ((lon + cp.pi) / (2 * cp.pi)) * pano.shape[1]
    v = ((cp.pi / 2 - lat) / cp.pi) * pano.shape[0]

    u = u.reshape(image_shape[0], image_shape[1]).astype(cp.float32)
    v = v.reshape(image_shape[0], image_shape[1]).astype(cp.float32)

    image = cp.empty((image_shape[0], image_shape[1], pano.shape[2]), dtype=pano.dtype)
    for i in range(pano.shape[2]):
        image[:, :, i] = ndimage.map_coordinates(pano[:, :, i], cp.stack([v, u]), order=1, mode='wrap')

    return image

问题根源分析 #

1. 旋转矩阵结构错误：
   原代码中R = cp.array([x, y, orientation])将基向量作为矩阵行存储，但标准3D旋转矩阵需以列向量构成正交矩阵，用于将局部相机坐标系的点转换到世界坐标系。行式矩阵会导致旋转方向完全颠倒，多分量朝向时误差被放大。
2. 矩阵乘法顺序错误：
   rotate_xyz中cp.matmul(R, xyz)的顺序对应行式矩阵左乘行向量，不符合坐标变换的数学定义，会导致点的旋转方向和轴出现偏差。
3. Up向量方向不匹配：
   默认up向量[0,-1,0]与等矩形全景图的坐标系（y轴正方向为天顶）相悖，会导致俯仰方向的转换出现翻转，进一步加剧畸变。

修正后的代码 #

import cupy as cp
from scipy import ndimage

def calculate_R(orientation, up=cp.array([0, 1, 0])):
    # 归一化朝向与up向量
    orientation = orientation / cp.linalg.norm(orientation)
    up = up / cp.linalg.norm(up)
    
    # 处理朝向与up平行的特殊情况
    if cp.allclose(orientation, up):
        return cp.eye(3)
    if cp.allclose(orientation, -up):
        return cp.array([[1, 0, 0], [0, -1, 0], [0, 0, -1]])
    
    # 构建正交基：x(右), y(上), z(朝向)
    z = orientation
    x = cp.cross(up, z)
    x = x / cp.linalg.norm(x)
    y = cp.cross(z, x)  # 确保y轴与z、x严格正交
    
    # 旋转矩阵：列向量构成，实现局部->世界坐标转换
    R = cp.array([x, y, z]).T
    
    return R

def rotate_xyz(R, x, y, z):
    # 将局部坐标点(N,3)与旋转矩阵相乘，转换为世界坐标
    xyz = cp.stack((x.ravel(), y.ravel(), z.ravel()), axis=1)
    xyz_world = cp.matmul(xyz, R)
    x, y, z = xyz_world[:,0].reshape(x.shape), xyz_world[:,1].reshape(y.shape), xyz_world[:,2].reshape(z.shape)
    return x, y, z

def transform_back(pano, orientation, image_shape, fov=[90, 90]):
    orientation = orientation / cp.linalg.norm(orientation)

    fov_x = cp.radians(fov[0])
    fov_y = cp.radians(fov[1])

    x_max = cp.tan(fov_x / 2)
    y_max = cp.tan(fov_y / 2)

    # 生成相机局部坐标系下的透视平面点
    x = cp.linspace(-x_max, x_max, num=image_shape[1], dtype=cp.float32)
    y = cp.linspace(-y_max, y_max, num=image_shape[0], dtype=cp.float32)
    x, y = cp.meshgrid(x, y)
    z = -cp.ones_like(x)  # 相机朝向为-z方向（看向场景）

    R = calculate_R(orientation)

    # 转换到世界坐标系
    x, y, z = rotate_xyz(R, x, y, z)
    xyz = cp.stack([x, y, z], axis=-1).reshape(-1, 3)

    # 计算全景图经纬度
    lon = cp.arctan2(xyz[:,0], xyz[:,2])
    lat = cp.arcsin(xyz[:,1] / cp.linalg.norm(xyz, axis=1))

    # 转换为全景图像素坐标
    u = ((lon + cp.pi) / (2 * cp.pi)) * pano.shape[1]
    v = ((cp.pi / 2 - lat) / cp.pi) * pano.shape[0]

    u = u.reshape(image_shape[0], image_shape[1]).astype(cp.float32)
    v = v.reshape(image_shape[0], image_shape[1]).astype(cp.float32)

    # 采样生成透视图像
    image = cp.empty((image_shape[0], image_shape[1], pano.shape[2]), dtype=pano.dtype)
    for i in range(pano.shape[2]):
        image[:, :, i] = ndimage.map_coordinates(pano[:, :, i], cp.stack([v, u]), order=1, mode='wrap')

    return image

关键修正说明 #

• 旋转矩阵重构：将基向量转为列向量构建矩阵，符合3D坐标变换的数学规范。
• 乘法顺序调整：改为点矩阵与旋转矩阵右乘，正确实现局部到世界坐标的转换。
• Up向量修正：改为[0,1,0]匹配全景图坐标系，避免俯仰方向翻转。
• 特殊情况优化：平行朝向时返回更合理的单位矩阵/翻转矩阵，减少错误分支。

内容的提问来源于stack exchange，提问作者Fred1313