使用GDAL库重新投影坐标

Reprojecting coordinates using GDAL library

本文关键字：投影 坐标 新投影 GDAL 使用      更新时间：2023-10-16

我正在尝试将WGS84系统中的lat/lon坐标重新投影到SIRGAS 2000坐标系中的UTM坐标。

使用GDAL，我开始将同一坐标系（29N）中的已知utm坐标更改为lat/lon对应的，只是为了检查我是否编写了正确的代码（此处省略错误检查）：

OGRSpatialReference monUtm;
monUtm.SetWellKnownGeogCS("WGS84");
monUtm.SetUTM(29, true);
OGRSpatialReference monGeo;
monGeo.SetWellKnownGeogCS("WGS84");
OGRCoordinateTransformation* coordTrans = OGRCreateCoordinateTransformation(&monUtm, &monGeo);
double x = 621921.3413490148;
double y = 4794536.070196861;
int reprojected = coordTrans->Transform(1, &x, &y);
// If OK, print the coords.
delete coordTrans;
coordTrans = OGRCreateCoordinateTransformation(&monGeo, &monUtm);
reprojected = coordTrans->Transform(1, &x, &y);
// If OK, Print the coords.
delete coordTrans;

坐标621921.3413490148, 4794536.070196861对应于加利西亚北部的蒙塞洛斯地区。前后转换似乎是正确的：lat/lon坐标是正确的，当投影回UTM时，我得到了与原始坐标相同的：

UTM:          621921.34135 , 4794536.0702
Lat/lon:      43.293779579 , -7.4970160261
Back to UTM:  621921.34135 , 4794536.0702

现在，从WGS84 lat/long重新投影到SIRGAS 2000 UTM:

// Rodovia dos Tamoios, Brazil:
//  - UTM          -> 23 S
//  - WGS 84       -> EPSG:32723
//  - SIRGAS 2000  -> EPSG:31983
OGRSpatialReference wgs;
wgs.SetWellKnownGeogCS("WGS84");
OGRSpatialReference sirgas;
sirgas.importFromEPSG(31983);
coordTrans = OGRCreateCoordinateTransformation(&wgs, &sirgas);
double x = -23.57014667;
double y = -45.49159617;
reprojected = coordTrans->Transform(1, &x, &y);
// If OK, print results
delete coordTrans;
coordTrans = OGRCreateCoordinateTransformation(&sirgas, &wgs);
reprojected = coordTrans->Transform(1, &x, &y);
// If OK, print results.

这并没有给出相同的结果：

WGS84 Lat/Lon input:      -23.57014667  ,  -45.49159617
SIRGAS 2000 UTM output:   2173024.0216  ,  4734004.2131
Back to WGS84 Lat/Lon:    -23.570633824 ,  -45.491627598

正如您所看到的，与第一个测试用例不同，最初的WGS84 lat/lon和back-to_WGS84 lat/lon坐标并不完全相同。此外，UTM x-coord有7位数字（我以为它被限制为6（？））。

在谷歌地图中，我们可以看到两个点之间有27米的差距（原始点用圆圈表示。我的"背部重投影"点用匕首代表）。

最后，问题是：我做的这个谴责是对的吗？如果是，为什么第二个测试用例中的重投影之间有27米的差异？

问题是，您需要将轴顺序交换为使用笛卡尔X/Y空间或Lon/Lat，而不是"Lat/Lon"顺序。

设置此项应该有效。

double x = -45.49159617;  // Lon
double y = -23.57014667;  // Lat

您从往返转换中看到的差异是由于交换轴顺序而投影到UTM区域的边界之外。

相关文章：