又是掉头发的一天，今天的任务是通过图片中心点的地理坐标以及图片中某点的像素坐标（即这个点位于图片中的x,y坐标）计算该点的地理/投影坐标。经过一整天的搜索，发现网上并没有这方面的教程。然后就是想啊想，头发一抓一大把，终于在网上零零散散的教程以及不断摸索中解决了这个问题。

大致思路就是，先获取图片相对真北方向的偏转角以及该点和图片中心的连线与图片的正北方向夹角；然后将图片中心点的地理坐标转换为投影坐标（如果这一步没有中心点的地理坐标，那么你就不用继续往下看了）；最后就是通过图片分辨率计算点到中心的实际距离，再通过夹角和中心点的投影坐标加加减减即可。话虽这么说，但实施起来真心不简单啊！！！

一、准备工作

1.获取图片中心点的经纬度坐标/投影坐标（必须）

如果只有一张图片的话，你完全可以右键图片，查看其属性，里面就有经纬度坐标。同时也可以使用Python实现，之前分享过【Python&GIS】判断图片中心点/经纬度点是否在某个面内，所以这里就不解释了，直接上代码。

def Get_LatLon(path_image):""":param path_image: 输入图片路径:return: 返回中心点经纬度"""# 获取图片的经纬度信息f = open(path_image, 'rb')contents = exifread.process_file(f)longitude = contents["GPS GPSLongitude"].valueslongitude_f = longitude[0].num/longitude[0].den + (longitude[1].num/longitude[1].den/60) + (longitude[2].num/longitude[2].den/3600)latitude = contents["GPS GPSLatitude"].valueslatitude_f = latitude[0].num/latitude[0].den + (latitude[1].num/latitude[1].den/60) + (latitude[2].num/latitude[2].den/3600)f.close()return longitude_f, latitude_f

2.获取图片与真北方向的偏转角（非必须）

如果你的图片并不是无人机拍的，或者拍摄时图片与真北方向并无夹角，那么就不需要这一步。这一步主要就是矫正相机拍摄时的偏转，每种无人机的参数名可能不一样，所以需要自行修改。

下面的代码也可以获取图片中心点的经纬度，但我之前以及用过第一步的代码了，所以就继续使用第一步的了，你们可以自己简化。

def Get_Image_Yaw_angle(file_path):""":param file_path: 输入图片路径:return: 图片的偏航角"""# 获取图片偏航角b = b"\x3c\x2f\x72\x64\x66\x3a\x44\x65\x73\x63\x72\x69\x70\x74\x69\x6f\x6e\x3e"a = b"\x3c\x72\x64\x66\x3a\x44\x65\x73\x63\x72\x69\x70\x74\x69\x6f\x6e\x20"img = open(file_path, 'rb')data = bytearray()dj_data_dict = {}flag = Falsefor line in img.readlines():if a in line:flag = Trueif flag:data += lineif b in line:breakif len(data) > 0:data = str(data.decode('ascii'))lines = list(filter(lambda x: 'drone-dji:' in x, data.split("\n")))for d in lines:d = d.strip()[10:]key, value = d.split("=")dj_data_dict[key] = value# print("Image_yaw",dj_data_dict["FlightYawDegree"][1:-1])return float(dj_data_dict["FlightYawDegree"][1:-1])

3.获取图片目标点与中心点的连线和图片正北方向的夹角（必须）

通俗来说，就是该点与图片正上方的夹角，同样这步也是用来矫正偏转的，获取该角度后，即可得到该点在投影坐标系中与真北方向的夹角。

def Yaw_angle_calculation(x, y, x_Center, y_Center):""":param x: 输入图片目标点的栅格坐标x:param y: 输入图片目标点的栅格坐标y:return: 目标点相对于中心点的偏转角度"""# 获取目标的偏转角if x == x_Center and y == y_Center:Yaw_angle = 0elif x == x_Center and y != y_Center:if y > y_Center:Yaw_angle = 180elif y < 1500:Yaw_angle = 0elif x != x_Center and y == y_Center:if x >x_Center:Yaw_angle = 90elif x < x_Center:Yaw_angle = 270else:if x > x_Center:if y < y_Center:Yaw_angle = math.degrees(math.atan((x - x_Center) / (y_Center - y)))else:Yaw_angle = 180 - math.degrees(math.atan((x - x_Center) / (y - y_Center)))else:if y < y_Center:Yaw_angle = 360 - math.degrees(math.atan((x_Center - x) / (y_Center - y)))else:Yaw_angle = 180 + math.degrees(math.atan((x_Center - x) / (y - y_Center)))# print("Yaw_angle",Yaw_angle)return Yaw_angle

4.将图片中心点的地理坐标转为投影坐标（必须）

一般来说图片的中心点都是GPS坐标，即WGS84地理坐标系的经纬度。但当我们计算距离时，需要使用到投影坐标系，所以这里需要将其转换一下。我这里是WGS84地理坐标系转为UTM 51N投影坐标系，你们视情况修改。

def LonLat_Meter(lat, lon):""":param lat: 输入WGS84经度:param lon: 输入WGS84纬度:return: 返回WGS84/UTM 51N的x,y"""source = osr.SpatialReference()# 初始化osr.SpatialReference对象以形成一个合法的坐标系统source.ImportFromEPSG(4326)# 向对象中写入Source_EPSG坐标系统target = osr.SpatialReference()target.ImportFromEPSG(32651)# 这里是用内置的EPSG对应的坐标系作为转换参数transform = osr.CoordinateTransformation(source, target)point = ogr.CreateGeometryFromWkt("POINT (%s %s)" % (lat, lon))# 报错的话，将经纬度翻过来point.Transform(transform)# print(point.GetX(), point.GetY())return point.GetX(), point.GetY()

二、转换函数

1.前期工作完成后，就可以实现转换了（必须！）

这里输入参数为图片偏转角、目标点的像素坐标（栅格坐标）、中心点的像素坐标、栅格分辨率（空间分辨率）、中心点的投影坐标。

def Target_Meter(Image_Yaw_angle, x, y, x_Center, y_Center, ratio,  meter_X, meter_Y):""":param Image_Yaw_angle: 输入图片的偏航角:param x: 输入图片目标点的栅格坐标x:param y: 输入图片目标点的栅格坐标y:param x_Center: 输入图片中心点的栅格坐标x:param y_Center: 输入图片中心点的栅格坐标y:param meter_X: 输入图片中心点的投影坐标x（UTM/WGS84 51N）:param meter_Y: 输入图片中心点的投影坐标y（UTM/WGS84 51N）:return: 目标点的投影坐标x,y"""if Image_Yaw_angle < 0:Image_Yaw_angle = 360 +Image_Yaw_angleYaw_angle = Yaw_angle_calculation(x, y, x_Center, y_Center)sum_angle = Image_Yaw_angle + Yaw_angleif sum_angle >= 360:sum_angle = sum_angle -360if sum_angle == 0:meter_x = meter_Xmeter_y = meter_Y - (y_Center-y)*ratioelif sum_angle == 90:meter_x = meter_X + (x-x_Center)*ratiometer_y = meter_Yelif sum_angle == 180:meter_x = meter_Xmeter_y = meter_Y + (y-y_Center)*ratioelif sum_angle == 270:meter_x = meter_X - (x_Center-x)*ratiometer_y = meter_Yelif sum_angle == 360:meter_x = meter_Xmeter_y = meter_Y - (x_Center-y)*ratioelif 0 < sum_angle < 90:meter_x = meter_X + math.sin(math.radians(sum_angle))*math.sqrt(math.pow(x-x_Center, 2)+math.pow(y-y_Center, 2))*ratiometer_y = meter_Y + math.cos(math.radians(sum_angle))*math.sqrt(math.pow(x-x_Center, 2)+math.pow(y-y_Center, 2))*ratioelif 90 < sum_angle < 180:meter_x = meter_X + math.sin(math.radians(180-sum_angle)) * math.sqrt(math.pow(x - x_Center, 2) + math.pow(y -y_Center, 2))*ratiometer_y = meter_Y - math.cos(math.radians(180-sum_angle)) * math.sqrt(math.pow(x - x_Center, 2) + math.pow(y -y_Center, 2))*ratioelif 180 < sum_angle < 270:meter_x = meter_X - math.sin(math.radians(sum_angle-180)) * math.sqrt(math.pow(x - x_Center, 2) + math.pow(y - y_Center, 2))*ratiometer_y = meter_Y - math.cos(math.radians(sum_angle-180)) * math.sqrt(math.pow(x - x_Center, 2) + math.pow(y - y_Center, 2))*ratioelif 270 < sum_angle < 360:meter_x = meter_X - math.sin(math.radians(360-sum_angle)) * math.sqrt(math.pow(x - x_Center, 2) + math.pow(y - y_Center, 2))*ratiometer_y = meter_Y + math.cos(math.radians(360-sum_angle)) * math.sqrt(math.pow(x - x_Center, 2) + math.pow(y - y_Center, 2))*ratioreturn meter_x, meter_y

2.将目标点的投影坐标转为地理坐标（非必须）

我这里因为项目最后需要的是经纬度坐标，所以在计算得到目标点的投影坐标后，还需要转换一下。你们视情况而定，如果你们只要投影坐标，那么第5步输出的就已经是投影坐标啦。

def Meter_LonLat(x, y):""":param x: 输入WGS84/UTM 51N的x:param y: 输入WGS84/UTM 51N的y:return: 返回WGS84的经纬度"""source = osr.SpatialReference()# 初始化osr.SpatialReference对象以形成一个合法的坐标系统source.ImportFromEPSG(32651)# 向对象中写入Source_EPSG坐标系统target = osr.SpatialReference()target.ImportFromEPSG(4326)# 这里是用内置的EPSG对应的坐标系作为转换参数transform = osr.CoordinateTransformation(source, target)point = ogr.CreateGeometryFromWkt("POINT (%s %s)" % (x, y))# 报错的话，将经纬度翻过来point.Transform(transform)# print(point.GetX(), point.GetY())return point.GetX(), point.GetY()

三、完整代码

# -*- coding: utf-8 -*-
"""
@Time ： 2023/1/5 17:45
@Auth ： RS迷途小书童
@File ：Picture coordinate system to geographic coordinate.py
@IDE ：PyCharm
@Purpose ：图片某点的像素坐标转为地理坐标系
"""
import math
import exifread
from osgeo import ogr, osrdef Get_LatLon(path_image):""":param path_image: 输入图片路径:return: 返回中心点经纬度"""# 获取图片的经纬度信息f = open(path_image, 'rb')contents = exifread.process_file(f)longitude = contents["GPS GPSLongitude"].valueslongitude_f = longitude[0].num/longitude[0].den + (longitude[1].num/longitude[1].den/60) + (longitude[2].num/longitude[2].den/3600)latitude = contents["GPS GPSLatitude"].valueslatitude_f = latitude[0].num/latitude[0].den + (latitude[1].num/latitude[1].den/60) + (latitude[2].num/latitude[2].den/3600)f.close()return longitude_f, latitude_fdef LonLat_Meter(lat, lon):""":param lat: 输入WGS84经度:param lon: 输入WGS84纬度:return: 返回WGS84/UTM 51N的x,y"""source = osr.SpatialReference()# 初始化osr.SpatialReference对象以形成一个合法的坐标系统source.ImportFromEPSG(4326)# 向对象中写入Source_EPSG坐标系统target = osr.SpatialReference()target.ImportFromEPSG(32651)# 这里是用内置的EPSG对应的坐标系作为转换参数transform = osr.CoordinateTransformation(source, target)point = ogr.CreateGeometryFromWkt("POINT (%s %s)" % (lat, lon))# 报错的话，将经纬度翻过来point.Transform(transform)# print(point.GetX(), point.GetY())return point.GetX(), point.GetY()def Meter_LonLat(x, y):""":param x: 输入WGS84/UTM 51N的x:param y: 输入WGS84/UTM 51N的y:return: 返回WGS84的经纬度"""source = osr.SpatialReference()# 初始化osr.SpatialReference对象以形成一个合法的坐标系统source.ImportFromEPSG(32651)# 向对象中写入Source_EPSG坐标系统target = osr.SpatialReference()target.ImportFromEPSG(4326)# 这里是用内置的EPSG对应的坐标系作为转换参数transform = osr.CoordinateTransformation(source, target)point = ogr.CreateGeometryFromWkt("POINT (%s %s)" % (x, y))# 报错的话，将经纬度翻过来point.Transform(transform)# print(point.GetX(), point.GetY())return point.GetX(), point.GetY()def Get_Image_Yaw_angle(file_path):""":param file_path: 输入图片路径:return: 图片的偏航角"""# 获取图片偏航角b = b"\x3c\x2f\x72\x64\x66\x3a\x44\x65\x73\x63\x72\x69\x70\x74\x69\x6f\x6e\x3e"a = b"\x3c\x72\x64\x66\x3a\x44\x65\x73\x63\x72\x69\x70\x74\x69\x6f\x6e\x20"img = open(file_path, 'rb')data = bytearray()dj_data_dict = {}flag = Falsefor line in img.readlines():if a in line:flag = Trueif flag:data += lineif b in line:breakif len(data) > 0:data = str(data.decode('ascii'))lines = list(filter(lambda x: 'drone-dji:' in x, data.split("\n")))for d in lines:d = d.strip()[10:]key, value = d.split("=")dj_data_dict[key] = value# print("Image_yaw",dj_data_dict["FlightYawDegree"][1:-1])return float(dj_data_dict["FlightYawDegree"][1:-1])def Yaw_angle_calculation(x, y, x_Center, y_Center):""":param x: 输入图片目标点的栅格坐标x:param y: 输入图片目标点的栅格坐标y:return: 目标点相对于中心点的偏转角度"""# 获取目标的偏转角if x == x_Center and y == y_Center:Yaw_angle = 0elif x == x_Center and y != y_Center:if y > y_Center:Yaw_angle = 180elif y < 1500:Yaw_angle = 0elif x != x_Center and y == y_Center:if x >x_Center:Yaw_angle = 90elif x < x_Center:Yaw_angle = 270else:if x > x_Center:if y < y_Center:Yaw_angle = math.degrees(math.atan((x - x_Center) / (y_Center - y)))else:Yaw_angle = 180 - math.degrees(math.atan((x - x_Center) / (y - y_Center)))else:if y < y_Center:Yaw_angle = 360 - math.degrees(math.atan((x_Center - x) / (y_Center - y)))else:Yaw_angle = 180 + math.degrees(math.atan((x_Center - x) / (y - y_Center)))# print("Yaw_angle",Yaw_angle)return Yaw_angledef Target_Meter(Image_Yaw_angle, x, y, x_Center, y_Center, ratio,  meter_X, meter_Y):""":param Image_Yaw_angle: 输入图片的偏航角:param x: 输入图片目标点的栅格坐标x:param y: 输入图片目标点的栅格坐标y:param x_Center: 输入图片中心点的栅格坐标x:param y_Center: 输入图片中心点的栅格坐标y:param meter_X: 输入图片中心点的投影坐标x（UTM/WGS84 51N）:param meter_Y: 输入图片中心点的投影坐标y（UTM/WGS84 51N）:return: 目标点的投影坐标x,y"""if Image_Yaw_angle < 0:Image_Yaw_angle = 360 +Image_Yaw_angleYaw_angle = Yaw_angle_calculation(x, y, x_Center, y_Center)sum_angle = Image_Yaw_angle + Yaw_angleif sum_angle >= 360:sum_angle = sum_angle -360if sum_angle == 0:meter_x = meter_Xmeter_y = meter_Y - (y_Center-y)*ratioelif sum_angle == 90:meter_x = meter_X + (x-x_Center)*ratiometer_y = meter_Yelif sum_angle == 180:meter_x = meter_Xmeter_y = meter_Y + (y-y_Center)*ratioelif sum_angle == 270:meter_x = meter_X - (x_Center-x)*ratiometer_y = meter_Yelif sum_angle == 360:meter_x = meter_Xmeter_y = meter_Y - (x_Center-y)*ratioelif 0 < sum_angle < 90:meter_x = meter_X + math.sin(math.radians(sum_angle))*math.sqrt(math.pow(x-x_Center, 2)+math.pow(y-y_Center, 2))*ratiometer_y = meter_Y + math.cos(math.radians(sum_angle))*math.sqrt(math.pow(x-x_Center, 2)+math.pow(y-y_Center, 2))*ratioelif 90 < sum_angle < 180:meter_x = meter_X + math.sin(math.radians(180-sum_angle)) * math.sqrt(math.pow(x - x_Center, 2) + math.pow(y -y_Center, 2))*ratiometer_y = meter_Y - math.cos(math.radians(180-sum_angle)) * math.sqrt(math.pow(x - x_Center, 2) + math.pow(y -y_Center, 2))*ratioelif 180 < sum_angle < 270:meter_x = meter_X - math.sin(math.radians(sum_angle-180)) * math.sqrt(math.pow(x - x_Center, 2) + math.pow(y - y_Center, 2))*ratiometer_y = meter_Y - math.cos(math.radians(sum_angle-180)) * math.sqrt(math.pow(x - x_Center, 2) + math.pow(y - y_Center, 2))*ratioelif 270 < sum_angle < 360:meter_x = meter_X - math.sin(math.radians(360-sum_angle)) * math.sqrt(math.pow(x - x_Center, 2) + math.pow(y - y_Center, 2))*ratiometer_y = meter_Y + math.cos(math.radians(360-sum_angle)) * math.sqrt(math.pow(x - x_Center, 2) + math.pow(y - y_Center, 2))*ratioreturn meter_x, meter_yif __name__ == "__main__":path = "G:/DJI_26_W.jpeg"x = 1600y = 2500x_Center = 2000y_Center = 1500ratio = 0.046longitude_f, latitude_f = Get_LatLon(path)meter_X, meter_Y = LonLat_Meter(latitude_f, longitude_f)Image_Yaw_angle = Get_Image_Yaw_angle(path)# 获取图片的偏转角meter_x, meter_y = Target_Meter(Image_Yaw_angle, x, y, x_Center, y_Center, ratio,  meter_X, meter_Y)lat_x, lon_y = Meter_LonLat(meter_x, meter_y)print(longitude_f, latitude_f)print(meter_X, meter_Y)print(Image_Yaw_angle)print(meter_x, meter_y)print(lat_x, lon_y)

大家可以自己添加for、while循环实现多张图片或图片中多个目标物的坐标转换。博主这里是为了目标识别后的目标物真实坐标的计算，所以将循环部分内置到了其他程序。

本文章主要是分享个人在学习Python过程中写过的一些代码。有些部分参考了前人以及官网的教程，如有侵权请联系作者删除，大家有问题可以随时留言交流，博主会及时回复。

【PythonGIS】无人机影像的像素坐标计算图片某点的地理/投影坐标相关推荐

图像坐标球面投影_地图开发知识之-投影坐标
地球投影由于地球是一个赤道略宽两极略扁的不规则的梨形球体,表面是一个不可展平的曲面,而地图通常是二维平面,因此在地图制图时首先要考虑把曲面转化成平面.然而,从几何意义上来说,球面是不可展平的曲面.要 ...
html热点区域确定坐标,html图片热点连接区域上的坐标是如何定位的？
通过和以及三个标签一起使用,可以在html页面中插入图片,并在该图片上创建一个或多个不同形状区域的热点链接,点击热点区域可以跳转到指定的其他页面.那么, 我们今天就学习一下,如何定位图片热点区域上的坐 ...
c++ 圆上任意点坐标计算_圆曲线上任意一点坐标计算
偏距 (m) 偏角 (° ) N(X) E(Y) DK07+590.000 3378604.933 453651.957 98°56′56.31″ 5 90 3378609.872 453652.73 ...
(Python)卫星RPC有理多项式模型读取与正反投影坐标计算原理与实现
(Python)卫星RPC有理多项式模型读取与正反投影坐标计算原理与实现文章目录 (Python)卫星RPC有理多项式模型读取与正反投影坐标计算原理与实现摘要 RPC几何定位模型介绍 RPC模型库 ...
点到直线 / 投影平面的坐标计算
以下为今天学习的笔记,内容包括:点到投影平面 / 直线的理论推导和代码部分 1. 点到投影平面 / 直线的理论推导以下内容为我的手写笔记.首先推导更直观的点到直线的投影坐标计算,然后电到投影平面的坐 ...
matlab附合导线坐标计算,“一步测量法”在农村土地确权测量中应用及精度分析...
"一步测量法"在农村土地确权测量中应用及精度分析作者:未知摘要:近几年来,随着国家城镇建设和新农村建设步伐加快,城镇地籍测量工作在全国范围内铺开.传统的平板测图已经很难满足工程 ...
python dataset[trans_科学网—Python GDAL 图像坐标，投影坐标，经纬度坐标三者映射及运行错误解决 - 吴妍潼的博文...
题记: 写该博客是因为自己经常遇到这个问题,而我发现网络上关于这方面浏览量高的一些代码竟然都有误,每次照搬都被虐得很惨.有一些同志在某些博客下方留言说代码有问题,而博主没有回应,也没有更改错误.为了自 ...
基于无人机影像结合Arcgis、Agisoft Metashape软件计算林地面积探讨
摘要:林地面积求算是森林资源调查.造林与采伐作业设计.生态系统环境损害(林业)司法鉴定等工作的重要工序,林地面积求算主要有仪器测量与图上量测两种方法,仪器测量方法主要有全站仪测算法.罗盘仪导线测量及G ...
基于Pix4Dmapper的运动结构恢复无人机影像三维模型重建
基于Pix4Dmapper的运动结构恢复无人机影像三维模型重建 1 背景知识 1.1 运动结构恢复方法原理 1.2 运动结构恢复方法流程 2 软件与数据准备 2.1 软件准备 2.2 数据准备 3 研 ...

【PythonGIS】无人机影像的像素坐标计算图片某点的地理/投影坐标