GDAL投影转换/坐标偏移笔记
1.OSR
# OSR(矢量数据投影) #作用:投影坐标系之间转换、地理坐标和投影坐标之间转换 #可用于几何对象和点(点属于几何对象) from osgeo import gdal import osr peters_sr = osr.SpatialReference() peters_sr.ImportFromProj4('...') ct = osr.CoordinateTransformation(web_mercator_sr, peters_sr) #world对象为web_mercator_sr投影,但没有分配srs world.Transform(ct) #对几何对象的转换 ct.TransformPoint(x, y) #对点的转换 #如果几何对象分配有srs,转换方法如下 world.TransformTo(web_mercator_sr) #读取投影的地理基准,用于和地理坐标进行转换 geosrs = peters_sr.CloneGeogCS()
2.PROJ
# PROJ(矢量数据投影) #作用:投影坐标之间转换、地理坐标和投影坐标之间转换 import pyproj utm_proj = pyproj.Proj('+proj=utm +zone=31 +ellps=WGS84') x, y = utm_proj(2.294, 48.858) #地理坐标转化为投影坐标 x1, y1 = utm_proj(x, y, inverse=True) #逆变换 wgs84 = pyproj.Proj('') nad27 = pyproj.Proj('') x, y = pyproj.transform(wgs84, nad27, 580744, 4504695) #投影坐标之间转换
3.栅格数据重投影
# 栅格数据重投影 # 栅格数据也可以重投影,但比矢量数据投影更复杂。 # 栅格数据需要处理栅格数据中像元会弯曲和移动的事实,一对一的映射并不存在 # 通常用最近邻域插值法、双线性插值和三次卷积插值法进行插值 # 方法:1.gdalwarp 2.AutoCreateWarpedVRT srs = osr.SpatialReference() srs.SetWellKnownGeogCS('WGS84') #UTM转无投影,即地理坐标系 old_ds = gdal.Open('nat_color.tif') vrt_ds = gdal.AutoCreateWarpedVRT(old_ds, None, srs.ExportToWkt(), gdal.GRA_Bilinear) #第一个默认值None,使用源栅格数据的srs;第二个如果None,表示不发生重投影 gdal.GetDriverByName('gtiff').CreateCopy('nat_color_wgs84.tif', vrt_ds) #该函数返回一个数据集对象,用CreateCopy函数保存为gtiff
4.GetGeoTransform
# GetGeoTransform,在真实坐标和栅格数据坐标具有相同srs情况下,计算坐标偏移。 #作用:图像坐标(行列号)和现实世界坐标(投影坐标或地理坐标)之间的转换。是仿射变换,不是投影转换,和上面不同。 #0、3 x\y坐标 起始点现实世界坐标 1、5 像素宽度和高度 2、4 x\y方向旋转角 gt = ds.GetGeoTransform() #正变换:图像坐标到现实世界坐标。正变换时输入行列号,输出的现实世界坐标是栅格图像srs下的坐标 inv_gt = gdal.InvGeoTransform(gt) #逆变换:现实世界坐标到图像坐标 offsets = gdal.ApplyGeoTransform(inv_gt, 465200, 5296000) #逆变换:输入的投影坐标具有和栅格图像的相同的srs xoff, yoff = map(int, offsets) #取整
5.gdal.Transformer
# gdal.Transformer,可计算相同srs下的坐标偏移;不能用于不同srs投影转换 #作用:现实世界坐标(投影坐标或地理坐标)与图像坐标(行列号)之间的转换、两个栅格图像之间像素坐标偏移(行列号),如镶嵌
#原理就是在相同srs情况下,计算图1的像素坐标到现实世界坐标的偏移,再从现实世界坐标偏移到图2的像素坐标。其实就是两次仿射变换(正、逆),从而把图1的像素坐标偏移到图2的像素坐标。
#所以不能用于不同srs情况,因为该函数没有内置不同srs的投影转换公式。只能用于相同srs下,两个栅格数据集坐标的偏移。 #这里in_ds和out_ds具有相同srs。转换目的是为了把不同栅格数据的图像坐标(行列号)进行偏移,方便镶嵌 trans = gdal.Transformer(in_ds, out_ds, []) #空白用于设置转换器选项 success, xyz = trans.TransformPoint(False, 0, 0) #False基于源数据计算目标栅格,反之为True。起始坐标为左上角 0,0 x,y,z = map(int, xyz) #图像坐标和现实世界坐标之间转换 trans = gdal.Transformer(out_ds, None, []) success, xyz = trans.TransformPoint(0, 1078, 648)
个人想法:
1.OSR \ PROJ \ AutoCreateWarpedVRT 分别用于矢量/点、栅格图像的投影转换,栅格与前两者有所不同。
2.GetGeoTransform是仿射变换,计算坐标偏移。(现实世界坐标(地理坐标或投影坐标)与图像坐标行列号),其中Geo[0/3]是现实世界坐标而不是行列号。
3.gdal.Transform本书中用于相同srs,不同栅格图像坐标(行列号)偏移的计算。书中用于镶嵌。
4.proj和osr都进行投影转换(可用于地理坐标和投影坐标之间的转换,是一种特例)。
5.计算栅格之间像素(行列号)偏移,这和仿射变换参数有关,而不是投影转换。
6.投影转换用OSR PROJ(矢量)/AutoCreateWarpedVRT(栅格);Geotransform用于相同srs,图像坐标和现实世界坐标转换;gdal.Transformer用于相同srs,栅格图像
间图像坐标偏移,也可以计算图像坐标和现实世界坐标转换。
7.osr.CoordinateTransformation是不同坐标系的投影转换;
而gdal.Transformer是相同srs下,不同栅格图像之间的图像坐标偏移(也可计算现实世界坐标和图像坐标之间的转换)。
两者可以组合使用,完成不同投影的栅格图像之间的图像坐标偏移。