GDAL之栅格重投影

来源:互联网 发布:ftp服务器软件 编辑:程序博客网 时间:2024/04/30 00:44

栅格重投影主要涉及到:空间坐标系的转换,栅格仿射变换系数的重新计算这两个主要方面。

示例代码是将一幅WGS84 UTM投影的影像重新投影到WGS84地理坐标系上。

#!/usr/bin/python# -*- coding: UTF-8 -*-import gdal,osrfrom gdalconst import *#获取源数据及栅格信息gdal.AllRegister()src_data = gdal.Open("smallaster.img")#获取源的坐标信息srcSRS_wkt=src_data.GetProjection()srcSRS=osr.SpatialReference()srcSRS.ImportFromWkt(srcSRS_wkt)#获取栅格尺寸src_width = src_data.RasterXSizesrc_height = src_data.RasterYSizesrc_count=src_data.RasterCount#获取源图像的仿射变换参数src_trans=src_data.GetGeoTransform()OriginLX_src=src_trans[0]OriginTY_src=src_trans[3]pixl_w_src=src_trans[1]pixl_h_src=src_trans[5]OriginRX_src=OriginLX_src+pixl_w_src*src_widthOriginBY_src=OriginTY_src+pixl_h_src*src_height#创建输出图像driver = gdal.GetDriverByName("GTiff")driver.Register()dst_data = driver.Create("tpix1.tif",src_width,src_height,src_count)#设置输出图像的坐标dstSRS=osr.SpatialReference()dstSRS.ImportFromEPSG(4326)#投影转换ct=osr.CoordinateTransformation(srcSRS,dstSRS)#计算目标影像的左上和右下坐标,即目标影像的仿射变换参数OriginLX_dst,OriginTY_dst,temp=ct.TransformPoint(OriginLX_src,OriginTY_src)OriginRX_dst,OriginBY_dst,temp=ct.TransformPoint(OriginRX_src,OriginBY_src)pixl_w_dst=(OriginRX_dst-OriginLX_dst)/src_widthpixl_h_dst=(OriginBY_dst-OriginTY_dst)/src_heightdst_trans=[OriginLX_dst,pixl_w_dst,0,OriginTY_dst,0,pixl_h_dst]#print outTransdstSRS_wkt=dstSRS.ExportToWkt()#设置仿射变换系数及投影dst_data.SetGeoTransform(dst_trans)dst_data.SetProjection(dstSRS_wkt)#重新投影gdal.ReprojectImage(src_data,dst_data,srcSRS_wkt,dstSRS_wkt,GRA_Bilinear)#创建四级金字塔gdal.SetConfigOption('HFA_USE_RRD', 'YES')dst_data.BuildOverviews(overviewlist=[2,4,8,16])#计算统计值for i in range(0,src_count):band=dst_data.GetRasterBand(i+1)band.SetNoDataValue(-99)print band.GetStatistics(0,1)


0 0
原创粉丝点击