GDAL setgeotransform()

Question

我一直在尝试使用GDAL库将地理信息添加到图像中。从GDAL网页的文档中，我可以发现我可以使用GDAlSetGeotransform()，为此，我需要一个六参数GDAL转换信息。在六个参数中，北上方向图像的x旋转和y旋转被认为是0。但在我的情况下，我没有北上图像。那么，如果我有图像的四个角坐标，我如何获得这些旋转值。或者，如果我有图像的四个角坐标，还有任何其他技术可以将地理信息添加到图像中。

Answer 1

我发现的最好的方法是使用最小二乘估计变换矩阵。它采用以下形式。

给定N (4+)组像素值及其匹配的世界坐标，构建一对变换。然后对变换进行求解。下面的示例以倍频程显示。如果使用C++，只需使用Eigen或OpenCV并构建摩尔-彭罗斯奇异值分解伪逆求解器。

#!/usr/bin/env octave -qf

%  Apply Geo Transform
function [Wx, Wy] = Apply_GeoTransform( Px, Py, adfGeoTransform )

    %  Multiply
    Wx = adfGeoTransform(1,2) * Px + adfGeoTransform(1,3) * Py +     adfGeoTransform(1,1);
    Wy = adfGeoTransform(1,5) * Px + adfGeoTransform(1,6) * Py + adfGeoTransform(1,4);

endfunction

%
%   Main Function
%   
%  Create Sample Inputs
W = [ 0, 50;
     50,  0;
      0,-50;
    -50,  0];

P = [ 0, 0;
     10, 0;
     10,10;
      0,10];

%  Build Input Transform (A)
Ai = [P(1,1), P(1,2), 1;
      P(2,1), P(2,2), 1;
      P(3,1), P(3,2), 1;
      P(4,1), P(4,2), 1];

%  Build Output Solution Pairs
Bx = [W(1,1);
      W(2,1);
      W(3,1);
      W(4,1)];

By = [W(1,2);
      W(2,2);
      W(3,2);
      W(4,2)];

%   Solve For Inverse x and y
Xi = Ai \ Bx;
Yi = Ai \ By;

%   Construct Actual Transform
adfGeoTransform = [Xi(3,1), Xi(1,1), Xi(2,1), Yi(3,1), Yi(1,1), Yi(2,1)];

%   Test Result
[Wx1, Wy1] = Apply_GeoTransform( 0, 0, adfGeoTransform );
[Wx2, Wy2] = Apply_GeoTransform(10, 0, adfGeoTransform );
[Wx3, Wy3] = Apply_GeoTransform(10,10, adfGeoTransform );
[Wx4, Wy4] = Apply_GeoTransform( 0,10, adfGeoTransform );
[Wx5, Wy5] = Apply_GeoTransform( 5, 5, adfGeoTransform );

printf('W1: %d, %d\n', Wx1, Wy1);
printf('W2: %d, %d\n', Wx2, Wy2);
printf('W3: %d, %d\n', Wx3, Wy3);
printf('W4: %d, %d\n', Wx4, Wy4);
printf('W5: %d, %d\n', Wx5, Wy5);