您现在的位置是:首页 >技术教程 >python将tif从wgs84转gcj02网站首页技术教程
python将tif从wgs84转gcj02
简介python将tif从wgs84转gcj02
# mark: 主要是dem从wgs84转成gcj02,步骤为wgs84的4326转成3857,之后投影进行偏移,再将偏移后的3857转成4326
from osgeo import gdal, osr
# 经纬度转投影
def tif4326To3857(input_file, output_file):
options = gdal.WarpOptions(format='GTiff', srcSRS='EPSG:4326', dstSRS='EPSG:3857', width=1452)
gdal.Warp(output_file, input_file, options=options)
print('转换完成')
# 投影转经纬度
def tif3857To4326(input_file, output_file):
options = gdal.WarpOptions(format='GTiff', srcSRS='EPSG:3857', dstSRS='EPSG:4326', width=1452)
gdal.Warp(output_file, input_file, options=options)
print('转换完成')
# wgs84投影偏移至gcj02
def tifWgs84ToGcj02(input_file, output_file):
dataset = gdal.Open(input_file)
im_width = dataset.RasterXSize #栅格矩阵的列数
im_height = dataset.RasterYSize # 栅格矩阵的列数
im_array = dataset.GetRasterBand(1).ReadAsArray()
im_projection = dataset.GetProjection()
im_GeoTransform = dataset.GetGeoTransform()
im_GeoTransform_new = (im_GeoTransform[0]+740.4591671, im_GeoTransform[1], im_GeoTransform[2], im_GeoTransform[3]+24.60901496, im_GeoTransform[4], im_GeoTransform[5])
# 创建新的tif文件
driver = gdal.GetDriverByName('GTiff')
out_tif = driver.Create(output_file, im_width, im_height, 1, gdal.GDT_Float32)
out_tif.SetGeoTransform(im_GeoTransform_new)
out_tif.SetProjection(im_projection) # 给新建图层赋予投影信息
# # 获取地理坐标系统信息,用于选取需要的地理坐标系统
# srs = osr.SpatialReference()
# # 导入tif文件的坐标系
# srs.ImportFromEPSG(3857) # 定义输出的坐标系为"WGS 84",AUTHORITY["EPSG","4326"]
# out_tif.SetProjection(srs.ExportToWkt()) # 给新建图层赋予投影信息
# 数据写出
out_tif.GetRasterBand(1).WriteArray(im_array) # 将数据写入内存,此时没有写入硬盘
out_tif.FlushCache() # 将数据写入硬盘
out_tif = None # 注意必须关闭tif文件
print('偏移完成')
if __name__ == '__main__':
file1 = 'C:\Users\ZJZN\Desktop\tiff\4326.tif'
file2 = 'C:\Users\ZJZN\Desktop\tiff\3857.tiff'
file3 = 'C:\Users\ZJZN\Desktop\tiff\3857-gcj02.tiff'
file4 = 'C:\Users\ZJZN\Desktop\tiff\4326-gcj02.tiff'
tif4326To3857(file1, file2)
tifWgs84ToGcj02(file2, file3)
tif3857To4326(file3, file4)
tif4326To3857(file1, file2) 是将wgs84的经纬度转成投影
tifWgs84ToGcj02(file2, file3) 是将投影的进行偏移
tif3857To4326(file3, file4) 再将偏移后的投影转成经纬度
注意:im_GeoTransform_new 中的偏移量是各地根据实际计算得到的。
gcj02-3857投影 | wgs84-3857投影 | 差值 | |
x | 12587949.13 | 12587208.67 | 740.4591671 |
y | 4337622.252 | 4337597.643 | 24.60901496 |
风语者!平时喜欢研究各种技术,目前在从事后端开发工作,热爱生活、热爱工作。