GDAL setgeotransform()

GDAL setgeotransform()

本文关键字:setgeotransform GDAL      更新时间:2023-10-16

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

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

给定N(4+)组像素值及其匹配的世界坐标,构建一对变换。然后求解变换。以下示例以八度音阶显示。如果使用C++,只需使用Eigen或OpenCV并构建Moore Penrose SVD伪逆解算器。

#!/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, %dn', Wx1, Wy1);
printf('W2: %d, %dn', Wx2, Wy2);
printf('W3: %d, %dn', Wx3, Wy3);
printf('W4: %d, %dn', Wx4, Wy4);
printf('W5: %d, %dn', Wx5, Wy5);