python+gdal地理坐标转投影坐标

1 前言

地理坐标系,是使用三维球面来定义地球表面位置,以实现通过经纬度对地球表面点位引用的坐标系。 地理坐标系经过地图投影操作后就变成了投影坐标系。而地图投影是按照一定的数学法则将地球椭球面上点的经维度坐标转换到平面上的直角坐标。

图片

2 流程

2.1 矢量数据地理坐标转投影坐标

要素feature的形状geometry是由一系列点坐标构成的,将每个要素的形状点一一进行坐标转换即可。

下面以WGS84坐标转UTM投影为例:

from osgeo import ogr,osr,gdal
import glob
import osdef vecter_WGS2UTM(shp_path, UTM_shp_path, longitude, is_north):ds = ogr.Open(shp_path)layer = ds.GetLayer(0)driver = ogr.GetDriverByName('ESRI Shapefile')# 创建输出文件if os.path.exists(UTM_shp_path):driver.DeleteDataSource(UTM_shp_path)out_ds = driver.CreateDataSource(UTM_shp_path)outlayer = out_ds.CreateLayer(UTM_shp_path[:-4],geom_type = ogr.wkbPolygon)# 当前地理参考spatialRef = layer.GetFeature(0).GetGeometryRef().GetSpatialReference()# 根据经度计算UTM区号,进而定义UTM投影zone = str(int(longitude/6) + 31)zone = int('326' + zone) if is_north else int('327' + zone)UTM_spatialRef = osr.SpatialReference()UTM_spatialRef.ImportFromEPSG(zone)# 投影转换coordinate_transfor = osr.CoordinateTransformation(spatialRef, UTM_spatialRef)# 定义输出属性表信息feature = layer.GetFeature(0)field_count = feature.GetFieldCount()field_names = []for attr in range(field_count):field_defn = feature.GetFieldDefnRef(attr)field_names.append(field_defn.GetName())outlayer.CreateField(field_defn)out_fielddefn = outlayer.GetLayerDefn()for feature in layer:geometry = feature.GetGeometryRef()geometry.Transform(coordinate_transfor)out_feature = ogr.Feature(out_fielddefn)out_feature.SetGeometry(geometry)for field_name in field_names:out_feature.SetField(field_name,feature.GetField(field_name))outlayer.CreateFeature(out_feature)feature.Destroy()out_feature.Destroy()# 清除缓存ds.Destroy()out_ds.Destroy()# 保存投影文件UTM_spatialRef.MorphFromESRI()prj_path = UTM_shp_path.replace(".shp",".prj")fn = open(prj_path,'w')fn.write(UTM_spatialRef.ExportToWkt())fn.close()

2.2 栅格数据地理坐标转投影坐标

栅格数据每个像素的地理/投影坐标是由仿射矩阵六参数和像素坐标计算得来的,所以先将仿射矩阵六参数进行转换,之后对栅格数据重采样即可。

下面以WGS84坐标转UTM投影为例:

def raster_WGS2UTM(raster_path, UTM_raster_path, longitude, is_north):raster_ds = gdal.Open(raster_path)raster_type = raster_ds.GetRasterBand(1).DataType# 栅格投影spatialRef = osr.SpatialReference()spatialRef.ImportFromWkt(raster_ds.GetProjection())# 根据经度计算UTM区号,进而定义UTM投影zone = str(int(longitude/6) + 31)zone = int('326' + zone) if is_north else int('327' + zone)UTM_spatialRef = osr.SpatialReference()UTM_spatialRef.ImportFromEPSG(zone)# 投影转换coordinate_transfor = osr.CoordinateTransformation(spatialRef, UTM_spatialRef)# 仿射矩阵六参数geotransform = raster_ds.GetGeoTransform()# 左上角upper left、右下角lower right坐标ul_x = geotransform[0]ul_y = geotransform[3]lr_x = geotransform[0]+geotransform[1]*raster_ds.RasterXSize+geotransform[2]*raster_ds.RasterYSizelr_y = geotransform[3]+geotransform[4]*raster_ds.RasterYSize+geotransform[5]*raster_ds.RasterYSize# 左上角、右下角在目标投影中的坐标(UTM_ul_x, UTM_ul_y, UTM_ul_z) = coordinate_transfor.TransformPoint(ul_y, ul_x)(UTM_lr_x, UTM_lr_y, UTM_lr_z) = coordinate_transfor.TransformPoint(lr_y, lr_x)# 创建目标图像文件driver = gdal.GetDriverByName("GTiff")UTM_raster_ds = driver.Create(UTM_raster_path, raster_ds.RasterXSize,raster_ds.RasterYSize,raster_ds.RasterCount,raster_type)# 转换后图像的分辨率resolution = (UTM_lr_x-UTM_ul_x)/raster_ds.RasterXSize# 转换后图像的六个放射变换参数UTM_transform = [UTM_ul_x, resolution, 0, UTM_ul_y, 0, -resolution]UTM_raster_ds.SetGeoTransform(UTM_transform)UTM_raster_ds.SetProjection(UTM_spatialRef.ExportToWkt())# 投影转换后需要做重采样gdal.ReprojectImage(raster_ds, UTM_raster_ds, spatialRef.ExportToWkt(), UTM_spatialRef.ExportToWkt(), gdal.GRA_Bilinear)# 关闭raster_ds = NoneUTM_raster_ds= None

来源:应用推广部

供稿:技术研发部

编辑:方梅

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处:http://www.mzph.cn/news/212935.shtml

如若内容造成侵权/违法违规/事实不符,请联系多彩编程网进行投诉反馈email:809451989@qq.com,一经查实,立即删除!

相关文章

基于STM32的四位数码管计数器设计与实现

✅作者简介:热爱科研的嵌入式开发者,修心和技术同步精进, 代码获取、问题探讨及文章转载可私信。 ☁ 愿你的生命中有够多的云翳,来造就一个美丽的黄昏。 🍎获取更多嵌入式资料可点击链接进群领取,谢谢支持!…

Docker Compose(容器编排)——9

目录 什么是 Docker Compose生活案例为什么要 Docker ComposeDocker Compose 的安装Docker Compose 的功能Docker Compose 使用场景Docker Compose 文件(docker-compose.yml) 文件语法版本文件基本结构及常见指令Docker Compose 命令清单 命令清单如下命…

垃圾回收器CMS和G1的区别

CMS和G1的区别 区别一: 使用范围不一样 CMS收集器是老年代的收集器,可以配合新生代的Serial和ParNew收集器一起使用 G1收集器收集范围是老年代和新生代。不需要结合其他收集器使用 区别二: STW的时间 CMS收集器以最小的停顿时间为目标的收…

C++11(下)

可变参数模板 C11的新特性可变参数模板能够创建可以接受可变参数的函数模板和类模板. 相比C98/03, 类模版和函数模版中只能含固定数量的模版参数, 可变模版参数无疑是一个巨大的改进, 然而由于可变模版参数比较抽象, 使用起来需要一定的技巧, 所以这块还是比较晦涩的.掌握一些基…

Vue 3项目的运行过程

概述: 使用Vite构建Vue 3项目后,当执行yarn dev命令启动服务时,项目就会运行起来,该项目会通过src\main.js文件将src\App.vue组件渲染到index.html文件的指定区域。 文件介绍: src\App.vue文件 Vue 3项目是由各种组件…

递归实现指数型枚举

title: 递归实现指数型枚举 date: 2023-12-10 19:29:20 tags: 递归 catgories: 算法进阶指南 —> 传送门 题目大意 从 1 ~ n n n 这 n n n 个整数随机选取任意多个,输出所有可能的选择方案 思路 这等价于每个整数可以选或者不选,所有的方案总数共有…

Spring Boot的日志

打印日志 打印日志的步骤: • 在程序中得到日志对象. • 使用日志对象输出要打印的内容 在程序中得到日志对象 在程序中获取日志对象需要使用日志工厂LoggerFactory,代码如下: package com.example.demo;import org.slf4j.Logger; import org.slf4j.LoggerFactory;public c…

STM32——继电器

继电器工作原理 单片机供电 VCC GND 接单片机, VCC 需要接 3.3V , 5V 不行! 最大负载电路交流 250V/10A ,直流 30V/10A 引脚 IN 接收到 低电平 时,开关闭合。

Go Fyne 入门

Fyne是一个用于创建原生应用程序的UI工具包,它简单易用,并且支持跨平台。以下是一个简单的Fyne教程,帮助你入门: 1. 安装Fyne 首先,确保你已经安装了Go语言。然后,在终端中运行以下命令来安装Fyne&#x…

android-xml语法

xml解析器 Android的XML文件语法是由Android系统中的解析器解析的。具体来说,Android使用了一个名为"Android Asset Packaging Tool (AAPT)"的工具来解析和处理XML文件。AAPT负责将XML文件编译为二进制格式,并在构建过程中将其打包到Android应…

第2节:Vue3 模板语法

Vue3 的模板语法主要包括以下几个部分&#xff1a; 插值表达式&#xff1a;使用双大括号 {{ }} 包裹变量&#xff0c;可以直接在模板中显示变量的值。 <div>{{ message }}</div>指令&#xff1a;以 v- 开头&#xff0c;后面跟一个自定义的名字&#xff0c;用于操…

从Centos-7升级到Centos-Stream-8

如果在正式环境升级&#xff0c;请做好数据备份以及重要配置备份&#xff01;因为升级会造一部分应用被卸载。 注意&#xff1a;升级前请备份好数据&#xff0c;升级可能会导致ssh的root用户无法登陆、网卡名称发生改变、引导丢失无法开机等问题。 1.安装epel源 yum -y install…

【Spring教程20】Spring框架实战:AOP(面对切面编程)知识总结

欢迎大家回到《Java教程之Spring30天快速入门》&#xff0c;本教程所有示例均基于Maven实现&#xff0c;如果您对Maven还很陌生&#xff0c;请移步本人的博文《如何在windows11下安装Maven并配置以及 IDEA配置Maven环境》&#xff0c;本文的上一篇为《利用 AOP通知获取数据代码…

软件测试(接口测试业务场景测试)

软件测试 手动测试 测试用例8大要素 编号用例名称&#xff08;标题&#xff09;模块优先级预制条件测试数据操作步骤预期结果 接口测试&#xff08;模拟http请求&#xff09; 接口用例设计 防止漏测方便分配工具&#xff0c;评估工作量和时间接口测试测试点 功能 单接口业…

华为OD机试真题-字符串变换最小字符串-2023年OD统一考试(C卷)

题目描述: 给定一个字符串s,最多只能进行一次变换,返回变换后能得到的最小字符串(按照字典序进行比较)。变换规则:交换字符串中任意两个不同位置的字符。 输入描述:一串小写字母组成的字符串s 输出描述:按照要求进行变换得到的最小字符串 补充说明:s是都是小写字符组成…

一台是阿里云,一台是腾讯云,一台是华为云,一台是百度云等多种公有云混合安装K8S集群

1. 修改主机名称和添加hosts #永久修改主机名 hostnamectl set-hostname master && bash #在master01上操作&#xff0c;阿里云服务器 hostnamectl set-hostname worker1 && bash #在node01上操作&#xff0c;阿里腾讯云服务器 hostnamectl set-ho…

利用Microsoft Visual Studio Installer Projects打包安装包

利用Microsoft Visual Studio Installer Projects打包安装包 具体步骤步骤1&#xff1a;安装扩展步骤2&#xff1a;创建 Setup 项目步骤3&#xff1a;设置属性步骤4&#xff1a;添加输出步骤5&#xff1a;添加文件步骤6&#xff1a;添加桌面快捷方式步骤7&#xff1a;添加菜单快…

【Table/SQL Api】Flink Table/SQL Api表转流读取MySQL

引入依赖 jdbc依赖 flink-connector-jdbc mysql-jdbc-driver 操作mysql数据库 <!-- Flink-Connector-Jdbc --><dependency><groupId>org.apache.flink</groupId><artifactId>flink-connector-jdbc_${scala.binary.version}</artifactId>…

Ubuntu上安装 Git

在 Ubuntu 上安装 Git 可以通过包管理器 apt 进行。以下是在 Ubuntu 上安装 Git 的步骤&#xff1a; 打开终端。你可以按 Ctrl Alt T 组合键来打开终端。 运行以下命令以确保你的系统的软件包列表是最新的&#xff1a; sudo apt update 安装 Git&#xff1a; sudo apt inst…

RT-DERT改进策略:AKConv即插即用,轻松涨点

摘要 提出了一种算法&#xff0c;用于生成任意尺寸卷积核的初始采样坐标。与常规卷积核相比&#xff0c;提出的AKConv实现了不规则卷积核的函数来提取特征&#xff0c;为各种变化目标提供具有任意采样形状和尺寸的卷积核&#xff0c;弥补了常规卷积的不足。在COCO2017和VisDro…