将utm坐标转换为参考相邻区域的坐标

Posted 2023-02-23

【中文标题】将utm坐标转换为参考相邻区域的坐标

【英文标题】：convert utm coordinate to coordinate in reference to neighboring zone

【发布时间】：2016-02-13 04:01:20

【问题描述】：

我在 NAD83 UTM 13N 中有数万个栅格。我正在尝试使用 arcpy.GetCellValue_management(raster.tif, point) 逐点提取数据，但数据的最西侧位于 UTM 12N 区域。有没有办法从 12N 获取坐标但参考 13N？项目要求是所有数据都在 UTM 13N 中，即使它是一个全州项目。我知道这很愚蠢。

【参考方案1】：

这可以通过 GDAL 来完成。将 dataPoints.shp 保存在所需的 UTM 网格中（即 UTM 13N），然后使用 GDAL 加载点图层、获取字段、获取几何、获取边界坐标、地理变换、栅格波段、点坐标（在 UTM 13N 中）和将栅格读取为数组。在所有栅格上建立一个循环，它的工作速度非常快。感谢Luke 提供详细信息here。

from osgeo import gdal, ogr
shp_filename = 'C:\\Path\\dataPoints_UTM13.shp'
ds = ogr.Open(shp_filename)
lyr = ds.GetLayer()
for feat in lyr:
    point_id_obj = feat.GetField("Sample")
    name = feat.GetField("Location_D")
    geom = feat.GetGeometryRef()
    mx, my = geom.GetX(), geom.GetY()
    path = 'C:\\RasterPath'
    raster = 'myraster'
    ras_open = gdal.Open('a\\b.tif'.format(a=path, b=raster))
    gt = aws_open.GetGeoTransform()
    rb = aws_open.GetRasterBand(1)
    px = abs(int((mx - gt[0]) / gt[1]))
    py = int((my - gt[3]) / gt[5])
    ras_obj = rb.ReadAsArray(px, py, 1, 1)
    print point_id_obj
    print name
    print mx, my