使用OpenCV、skimage与SciPy执行仿射变换结果不一致问题排查

场景描述 #

需要将单应性矩阵应用于50kx50k像素的超大图像，计划通过SciPy结合Dask实现高效处理，但小尺度测试中发现SciPy的变换结果与OpenCV、Skimage不符：SciPy能正确旋转图像，但缩放效果错误；反转矩阵后缩放正常，但旋转方向相反。

核心原因 #

不同库对变换方向、坐标系统定义的处理逻辑存在差异：

1. 变换方向差异：
   • OpenCV warpAffine 是正向映射（将原图像素映射到目标图）；Skimage transform.warp 是反向映射（为目标图像素找原图对应点，因此需要传入矩阵的逆）。
   • SciPy ndimage.affine_transform 默认采用反向映射，但矩阵作用的是数组的索引坐标（行优先，即第一个维度是y轴，第二个是x轴），而非图像的常规(x,y)坐标（x轴向右，y轴向下，左上角为原点）。
2. 矩阵维度适配问题：
   图像的(x,y)坐标对应数组的(y,x)索引，因此需要对仿射矩阵的线性变换部分进行转置+维度交换，平移量也需交换x、y顺序，才能匹配SciPy的索引空间。

修正后的代码示例 #

import scipy
import cv2
import skimage
import numpy as np
import matplotlib.pyplot as plt

# 创建变换矩阵
translation = np.array([[1, 0, 40],  # x方向平移40，y方向平移-48
                        [0, 1, -48],
                        [0, 0, 1]])

rotation_angle = np.deg2rad(45)  
rotation = np.array([[np.cos(rotation_angle), -np.sin(rotation_angle), 0],
                     [np.sin(rotation_angle), np.cos(rotation_angle), 0],  # 修正原代码错误：此处平移量应为0
                     [0, 0, 1]])

scaling = np.array([[2, 0, 0],  # x、y方向均缩放2倍
                    [0, 2, 0],
                    [0, 0, 1]])

# 组合变换（矩阵乘法顺序：缩放→旋转→平移）
affine_matrix = translation @ rotation @ scaling

# 加载测试图像
image = skimage.data.astronaut()
h, w = image.shape[:2]

# Skimage 变换结果
skimage_result = skimage.transform.warp(
        image,
        np.linalg.inv(affine_matrix),
        output_shape=(h, w)
    )

# OpenCV 变换结果
cv2_result = cv2.warpAffine(image, 
                        affine_matrix[:2],
                        dsize=(w, h)
                        )

# 修正SciPy变换调用
# 1. 反向映射需使用矩阵的逆
inv_affine = np.linalg.inv(affine_matrix)
# 2. 转换矩阵适配SciPy的索引坐标（图像(x,y) → 数组(y,x)）
linear_transform = inv_affine[:2, :2][[1,0], :][:, [1,0]]  # 交换行列+转置
offset = inv_affine[:2, 2][[1,0]]  # 交换x、y平移量

scipy_result_corrected = []
for ch in range(image.shape[2]):
    result = scipy.ndimage.affine_transform(
            image[..., ch],
            linear_transform,
            offset=offset,
            output_shape=(h, w),
            order=0,
        )
    scipy_result_corrected.append(result)
scipy_result_corrected = np.stack(scipy_result_corrected).transpose(1, 2, 0)

# 可视化对比
fig, ax = plt.subplots(1, 4, figsize=(16, 8))
ax[0].imshow(image)
ax[0].set_title("原图")
ax[1].imshow(skimage_result)
ax[1].set_title("skimage.transform.warp")
ax[2].imshow(cv2_result)
ax[2].set_title("cv2.warpAffine")
ax[3].imshow(scipy_result_corrected)
ax[3].set_title("修正后的SciPy变换结果")

plt.tight_layout()
plt.show()

关键修正点 #

1. 修正旋转矩阵错误：原代码中旋转矩阵的第二行第三列误写为1，改为0，消除额外平移干扰。
2. 适配SciPy坐标系统：将图像(x,y)坐标对应的矩阵转换为数组(y,x)索引对应的矩阵，通过交换行列+转置实现。
3. 明确反向映射逻辑：SciPy与Skimage一致，需使用变换矩阵的逆，同时指定output_shape保证输出尺寸与原图一致。

修正后，SciPy的变换结果将与OpenCV、Skimage完全对齐，可直接结合Dask处理超大尺寸图像。

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