用 Python 进行仿射变换答案

【问题标题】：Affine transformation of with Python用 Python 进行仿射变换
【发布时间】：2022-01-17 14:58:27
【问题描述】：

从 MATLAB 翻译我的代码我正在尝试使用 Python 3 应用仿射变换来对 TERRA-ASTER 卫星图像进行地理配准，

import numpy as np
import matplotlib.pyplot as plt
from skimage import transform as tf
from skimage import io, exposure

naivasha_rgb = io.imread('naivasha.tif')

使用来自四个场景角的坐标对 (src,dst) 的变换矩阵。

# upper left
# lower left
# upper right
# lower right
movingpoints = np.array([
    [36.214332,-0.319922],
    [36.096003,-0.878267],
    [36.770406,-0.400443],
    [36.652213,-0.958743]])
fixedpoints = np.array([
    [0,0],
    [0,4200],
    [4100,0],
    [4100,4200]])

使用 Skimage 包中的函数

tform = tf.estimate_transform('affine',movingpoints,fixedpoints)
newnaivasha_rgb = tf.warp(naivasha_rgb,tform.inverse)

tform 与 MATLAB 中的相同，但经过转置，但变形会创建平面绿色图像，而不是 MATLAB 中的预期图像（请参阅下面的 WeTransfer 链接）。

newnaivasha_rgb = exposure.rescale_intensity(newnaivasha_rgb,
    out_range='uint8')
    
plt.figure()
io.imshow(newnaivasha_rgb)
plt.axis('off')
plt.show()

有什么想法吗？根据我搜索的解决方案，数值表示低对比度或没有图像，具体取决于我使用 tform 还是 tform.inverse 和 tf.warp。

从WeTransfer (50 MB) 下载图像，包括输入图像 naivasha.tif 和输出图像 naivasha_georef_python.png 和 naivasha_georef_matlab.jpg。

【问题讨论】：

并不是我对skimage 包有很多了解，但我想知道您的moving_points 是否需要与naivasha_rgb 的np.array 索引相关而不是纬度/经度.
可能，至少我得到了一张图片——谢谢！所以 Python 函数的工作方式似乎与我在 MATLAB 中使用的不同。
我不知道skimage 是否会使用它们，但那个 TIFF 似乎没有任何 GeoTIFF 元数据。
没错，但 TIFF 不是 GeoTIFF。它是从原始的 TERRA-ASTER HDF 文件中处理和保存的，与简单的 MATLAB 解决方案相比，这在使用 Python 时已经是一个相当大的挑战。 mres.uni-potsdam.de/index.php/2017/02/23/…
我只提到了 GeoTIFF，因为 Python 无法从 naivasha.tif 确定纬度/经度。如果 Matlab 更方便，为什么不直接使用呢？

标签： satellite terra aster

【解决方案1】：

感谢 rickhg12hs 的建议，我现在能够解决问题。我将经度和纬度值转换为相对于图像的像素坐标，然后仿射变换产生预期的结果。谢谢 rickhg12hs！

import numpy as np
import matplotlib.pyplot as plt
from skimage import transform as tf
from skimage import io

#%%
naivasha_rgb = io.imread('naivasha.tif')

#%%
plt.figure()
io.imshow(naivasha_rgb)
plt.axis('off')
plt.show()

#%%
UL = [-0.319922,36.214332]
LL = [-0.878267,36.096003]
UR = [-0.400443,36.770406]
LR = [-0.958743,36.652213]

#%%
# 4200 rows correspond to the latitudes  [1]
# 4100 cols correspond to the longitudes [0]
movingpoints = np.array([
    [0.,4100.*(UL[1]-LL[1])/(UR[1]-LL[1])],
    [4200.*(UL[0]-LL[0])/(UL[0]-LR[0]),0.],
    [4200.*(UL[0]-UR[0])/(UL[0]-LR[0]),4100.],
    [4200.,4100.*(LR[1]-LL[1])/(UR[1]-LL[1])]])

#%%
fixedpoints = np.array([
    [0.,0.],
    [4200.,0.],
    [0.,4100.],
    [4200.,4100.]]) 

#%% 
fixedpoints = np.fliplr(fixedpoints)
movingpoints = np.fliplr(movingpoints)
    
#%%
tform = tf.estimate_transform('affine',movingpoints,fixedpoints)

#%%
newnaivasha_rgb = tf.warp(naivasha_rgb,tform)

#%%
plt.figure()
plt.imshow(np.flipud(newnaivasha_rgb),
    extent=[36.096003,36.770406,-0.319922,-0.958743])
ax=plt.gca()
ax.invert_yaxis()
plt.savefig('newnaivasha_rgb.png')
plt.show()

【讨论】：