博主准研究僧一枚，假期在老师指导下接触项目。
本博文可作为坐标转换，特别是布尔莎七参数法的学习资料。其python源码注释充分，也可作为python的学习项目。
程序UI界面如下，由于是自用程序，博主对美化UI不感兴趣，ui部分源码注释充分，颜控可自行修改调整。
PS:克总信徒可以私信博主，富坦！
目前的大地测量实践中，由于GSC2000坐标系的推广，需要做大量将54坐标与80坐标转换为2000坐标的工作。实际工作中，一般通过mapgis实现坐标转换。但是，在某些涉密部门，仍然需要编制自用的程序进行坐标转换。
博主查阅资料，发现这方面资料要么比较老，时效性差，要么泛泛而谈，特别是介绍转换的详细算法的很少。所以，博主写下这篇短文，简要介绍下布尔莎模型及最小二乘思想在参数反解中的应用。同时博主注意到尚无利用python实现布尔莎七参数转换的博主，亦将展示源码，源码注释充分，以供交流学习。
坐标转换中用到的正是布尔莎-沃尔夫（Bursa-Wolf）模型，一般用七参数法进行坐标转换。涉及到的七个参数为：三个位移参数，三个旋转参数，一个尺度参数。不知为何写不了公式，放图：
可转换为：
七参数的矩阵形式即为：
上面三个是位移参数，然后是三个旋转参数，最后k为尺度参数即缩放因子。
我们定义伪观测数据V：
在大多数情况下我们是不知道七个参数的具体值的，所以一般通过最小二乘法拟合反解出七个参数。按照最小二乘原理，实际数据即是伪观测数据与线性公式得到的数据之间误差的平方和最小时的参数。一般需要至少三个公共点的坐标值反解七参数。以三个公共点为例，定义：
使Q最小：
附上python代码如下：
Java留下的习惯，简单的主程序代码：
main.py
from ui import Application
from transform import Coordinate
app=Application()
app.master.title("基于布尔莎模型的坐标转换工具")#窗口标题
#编程过程中间检测性质的代码段
# test=Coordinate([1,2,3],[1,2,3],[1,2,3],[1,2,3],[1,2,3],[1,2,3])
# print(test.V)
# print(test.K)
# print(test.K.shape)
#主消息循环
app.mainloop()
其实写一个函数模块逻辑上更清晰，但出于学习的作用，写了一个坐标系类，带有可以自己转换自己的方法（大雾.......）
numpy真的是个很好用的模块！从某种程度上numpy已经取代了matlab在我干活中的作用。
transform.py
#用python利用最小二乘原理计算矩阵，就得引入numpy模块
from numpy import *
class Coordinate():
"""坐标Coordinate类"""
def __init__(self,X1,Y1,Z1,X2,Y2,Z2,sign,n):
"""规定一些坐标类的属性"""
if sign==False:
V1=mat([[X1[0]-X2[0]],[Y1[0]-Y2[0]],[Z1[0]-Z2[0]]])
V2=mat([[X1[1]-X2[1]],[Y1[1]-Y2[1]],[Z1[1]-Z2[1]]])
V3=mat([[X1[2]-X2[2]],[Y1[2]-Y2[2]],[Z1[2]-Z2[2]]])
#创建伪观测值矩阵V
self.V=vstack((V1,V2,V3))#按列合并
#hstack按行合并
#创建K'''矩阵
K1=mat([[1,0,0,0,-Z1[0],Y1[0],X1[0]],[0,1,0,Z1[0],0,-X1[0],Y1[0]],[0,0,1,-Y1[0],X1[0],0,Z1[0]]])
K2=mat([[1,0,0,0,-Z1[1],Y1[1],X1[1]],[0,1,0,Z1[1],0,-X1[1],Y1[1]],[0,0,1,-Y1[1],X1[1],0,Z1[1]]])
K3=mat([[1,0,0,0,-Z1[2],Y1[2],X1[2]],[0,1,0,Z1[2],0,-X1[2],Y1[2]],[0,0,1,-Y1[2],X1[2],0,Z1[2]]])
self.K=vstack((K1,K2,K3))
#可以用shape方法（e.shape）查看矩阵维数，以确定代码正确
#根据输入的公共点进行最小二乘反解
else:
k=1
self.V=mat([[X1[0]-X2[0]],[Y1[0]-Y2[0]],[Z1[0]-Z2[0]]])
while k<n:
V1=mat([[X1[k]-X2[k]],[Y1[k]-Y2[k]],[Z1[k]-Z2[k]]])
self.V=vstack((self.V,V1))
k+=1
k=1
self.K=mat([[1,0,0,0,-Z1[0],Y1[0],X1[0]],[0,1,0,Z1[0],0,-X1[0],Y1[0]],[0,0,1,-Y1[0],X1[0],0,Z1[0]]])
while k<n:
K1 = mat([[1, 0, 0, 0, -Z1[k], Y1[k], X1[k]], [0, 1, 0, Z1[k], 0, -X1[k], Y1[k]],
[0, 0, 1, -Y1[k], X1[k], 0, Z1[k]]])
self.K=vstack((self.K,K1))
k+=1
self.X=((self.K.T*self.K).I)*(self.K.T)*self.V#X矩阵即是布尔莎模型根据最小二乘原理得到的转换参数矩阵
#X矩阵的排列方式是三个位移参数三个旋转参数和最后一个缩放系数
def transform(self,x,y,z):
""""类自带的方法"""
Xa=mat([x,y,z])
Xa=Xa.T
K=mat([[1,0,0,0,-z,y,x],[0,1,0,z,0,-x,y],[0,0,1,-y,x,0,z]])
#Xb作为类的一个属性，调用transform方法，得到该属性
self.Xb=Xa+K*self.X
最后是程序的主干，包含ui部分，主要采用tkinter模块，excel的读写上采用xlrd/xlwt等模块，可做新手学习的实例：
ui.py
from transform import Coordinate
#导入tkinter模块
from tkinter import *
import tkinter.messagebox as messagebox
import tkinter.filedialog
#sys模块对程序进行基本控制
import sys
#excel读写两个模块
import xlrd
import xlwt
class Application(Frame):#创建一个容器，为所有Widget的父容器，继承Frame类
def __init__(self,master=None):#master是窗口管理器，用于管理窗口部件，比如按钮标签等等等，master为None即是自己管理自己，为顶级窗口
"""创建主窗口也是父容器"""
Frame.__init__(self,master)
self.pack()#这里图省事，主窗口采用pack布局，其实个人觉得grid布局更实用
self.creatWidgets()
def creatWidgets(self):
"""主程序"""
#程序标签语，采用grid布局，覆盖两列
self.firstLabel=Label(self,text='基于布尔莎模型的坐标转换程序'+'\n'+'请先输入三个公共点的坐标或选择公共点excel表进行拟合')
self.firstLabel.grid(row=0,columnspan=10)
#模式转换按钮用checkbutton控件来写，学习下checkbutton的用法
self.b=tkinter.BooleanVar()
#self.b是该类的扩展的重要属性，用以选择是还按之前的程序目标通过手动输入三个公共点坐标来拟合，还是通过读取公共点的excel表，读取n个公共点
#同时扩展属性self.n作为一个可变变量，即公共点个数
self.n=tkinter.IntVar()
self.modelButton=Checkbutton(self,variable=self.b,text='excel输入公共点')
self.modelButton.grid(row=7,column=8)
#输入框行注释
self.rowoneLabel=Label(self,text='转换前的公共测量点坐标')
self.rowoneLabel.grid(row=1,column=1,columnspan=3)
self.rowtwoLabel=Label(self,text='转换后的公共测量点坐标')
self.rowtwoLabel.grid(row=1,column=6,columnspan=3)
#输入框列注释
self.coloneLabel=Label(self,text='X: ')
self.coloneLabel.grid(row=2,column=0)
self.coltwoLabel=Label(self,text='Y: ')
self.coltwoLabel.grid(row=3,column=0)
self.colthreeLabel=Label(self,text='Z: ')
self.colthreeLabel.grid(row=4,column=0)
#创建旧坐标输入框
self.e11=Entry(self,width=15)#Xa
self.e11.grid(row=2,column=1)
self.e21=Entry(self,width=15)#Ya
self.e21.grid(row=3,column=1)
self.e31=Entry(self,width=15)#Za
self.e31.grid(row=4,column=1)
self.e12=Entry(self,width=15)#Xa
self.e12.grid(row=2,column=2)
self.e22=Entry(self,width=15)#Ya
self.e22.grid(row=3,column=2)
self.e32=Entry(self,width=15)#Za
self.e32.grid(row=4,column=2)
self.e13 = Entry(self,width=15) # Xa
self.e13.grid(row=2, column=3)
self.e23 = Entry(self,width=15) # Ya
self.e23.grid(row=3, column=3)
self.e33 = Entry(self,width=15) # Za
self.e33.grid(row=4, column=3)
#创建新坐标输入框
self.n11=Entry(self,width=15)#Xa
self.n11.grid(row=2,column=6)
self.n21=Entry(self,width=15)#Ya
self.n21.grid(row=3,column=6)
self.n31=Entry(self,width=15)#Za
self.n31.grid(row=4,column=6)
self.n12=Entry(self,width=15)#Xa
self.n12.grid(row=2,column=7)
self.n22=Entry(self,width=15)#Ya
self.n22.grid(row=3,column=7)
self.n32=Entry(self,width=15)#Za
self.n32.grid(row=4,column=7)
self.n13 = Entry(self,width=15) # Xa
self.n13.grid(row=2, column=8)
self.n23 = Entry(self,width=15) # Ya
self.n23.grid(row=3, column=8)
self.n33 = Entry(self,width=15) # Za
self.n33.grid(row=4, column=8)
#创建一块空白部分作为分割
self.middleLabel=Label(self,text=' ',width=15)
self.middleLabel.grid(row=2,column=4,columnspan=2,rowspan=3)
#给显示窗口添加标签
self.delayLabel=Label(self,text='布尔莎模型七参数矩阵拟合结果[Tx,Ty,Tz,wx,wy,wz,k]')
self.delayLabel.grid(row=5,columnspan=10)
#创建显示窗口，显示七参数的拟合结果
self.Xresult=StringVar()#利用StrtingVat储存值
self.delay1=Entry(self,width=120,state='readonly',text=self.Xresult)
self.delay1.grid(row=6,columnspan=10)
#创建拟合按钮
self.XButton=Button(self,text='拟合',command=self.first_step,width=15)
self.XButton.grid(row=7,column=5)
#创建三个路径显示的StringVar
self.file_name_read=StringVar()
self.file_name_write=StringVar()
self.public_point=StringVar()
#创建四个按钮
self.openButton=Button(self,text='打开',width=15,command=self.read_location)
self.openButton.grid(row=8,column=0)
self.writeButton=Button(self,text='路径',width=15,command=self.write_location)
self.writeButton.grid(row=9,column=0)
self.openButton1=Button(self,text='打开公共点excel表',width=15,command=self.read_public_point)
self.openButton1.grid(row=10,column=0)
self.creatButton=Button(self,text='输出',command=self.write_excel,width=15)
self.creatButton.grid(row=11,column=5)
#路径窗口
self.delay2=Entry(self,width=100,state='readonly',text=self.file_name_read)
self.delay2.grid(row=8,column=1,columnspan=8)
self.delay3=Entry(self,width=100,state='readonly',text=self.file_name_write)
self.delay3.grid(row=9,column=1,columnspan=8)
self.delay4=Entry(self,width=100,state='readonly',text=self.public_point)
self.delay4.grid(row=10,column=1,columnspan=8)
#作者介绍
self.writerLabel=Label(self,text='制作人——田睿，密斯卡托尼克大学助教，阿卡姆镇巡更人，印斯茅斯荣誉市民，拉莱耶城的监视者'+'\n'+'欢迎克苏鲁神话爱好者私信骚扰')
self.writerLabel.grid(row=12,columnspan=10,rowspan=2)
def first_step(self):
"""该函数主要进行第一步的七参数拟合工作"""
#为Coordinate类创建输入矩阵量
#由于经常会发生输入异常，这里写了一个异常处理代码
if self.b.get()==False:
# 由于经常会发生输入异常，这里写了一个异常处理代码
try:
X1=[float(self.e11.get()),float(self.e12.get()),float(self.e13.get())]
Y1=[float(self.e21.get()),float(self.e22.get()),float(self.e23.get())]
Z1=[float(self.e31.get()),float(self.e32.get()),float(self.e33.get())]
X2=[float(self.n11.get()),float(self.n12.get()),float(self.n13.get())]
Y2=[float(self.n21.get()),float(self.n22.get()),float(self.n23.get())]
Z2=[float(self.n31.get()),float(self.n32.get()),float(self.n33.get())]
#得到转换参数结果矩阵，用Xresult接收（已经定义为StringVar）
#cd 作为Application类的一个属性，接收Coordinate类的实例，这样写是方便后边的方法调用
#这样写不会存在风险，因为正常使用程序，必须先调用first_step方法
self.cd=Coordinate(X1,Y1,Z1,X2,Y2,Z2,self.b.get(),self.n.get())
X1=str(self.cd.X)#得到矩阵的字符串形式，但需要进行加工
#对字符串进行加工
X1=X1[2:-3]
X1=X1.replace(']',',')
X1=X1.replace('[', '')
X2='['+X1+']'
self.Xresult.set(X2)
except ValueError:
messagebox.showinfo('错误','如果要使用输入的三个公共点进行拟合的话，请将三个公共测量点的坐标输全，再拟合反解参数！')
else:
messagebox.showinfo('提示','拟合完毕，可以进行下一步工作')
else:
#！！！！！！！！！！！！！！！！！！
#这里强调一个错误，之前这样写的：read_public_excel()，出错NameError: name 'read_public_excel' is not defined
#在非结构方法中（java的叫法）引用必须用self. 记得java是this？？？好久不用Java忘了
self.read_public_excel()
def read_public_point(self):
file_name=tkinter.filedialog.askopenfilename()
self.public_point.set(file_name)
def read_location(self):
file_name=tkinter.filedialog.askopenfilename()
self.file_name_read.set(file_name)
def write_location(self):
file_name=tkinter.filedialog.askopenfilename()
self.file_name_write.set(file_name[:-1])
def read_public_excel(self):
"""读取公共点excel表"""
file_name=self.public_point.get()
data=xlrd.open_workbook(file_name)
table=data.sheets()[0]
self.n.set(table.nrows)#n为行数,即公共点个数
n1=self.n.get()
if n1<3:
messagebox.showinfo('错误', '拟合需要至少三个以上公共点的数据')
k=0
#定义输入的列表
X1=[]
Y1=[]
Z1=[]
X2=[]
Y2=[]
Z2=[]
while k<n1:
X1.append(float(table.cell_value(k,0)))
Y1.append(float(table.cell_value(k, 1)))
Z1.append(float(table.cell_value(k, 2)))
X2.append(float(table.cell_value(k,3)))
Y2.append(float(table.cell_value(k, 4)))
Z2.append(float(table.cell_value(k, 5)))
k+=1
cd=Coordinate(X1,Y1,Z1,X2,Y2,Z2,True,n1)
X1 = str(cd.X) # 得到矩阵的字符串形式，但需要进行加工
# 对字符串进行加工
X1 = X1[2:-3]
X1 = X1.replace(']', ',')
X1 = X1.replace('[', '')
X2 = '[' + X1 + ']'
self.Xresult.set(X2)
def read_excel(self):
"""读取excel表函数"""
file_name=self.file_name_read.get() #想获取stringvar中的内容为字符串，用get方法，用set存储
#读取excel
# data=xlrd.open_workbook('E:/workspace/BJ54transform_project/Data/test.xlsx')
data=xlrd.open_workbook(file_name)
table=data.sheets()[0]#获取其中一个工作表
i=table.nrows#i为行数
k=0#计数变量
#采用一个字典储存值
datas={'x':0,'y':0,'z':0}
saves=[]#存储字典单元的列表
while k<i:
datas['x']=table.cell(k,0).value
datas['y']=table.cell(k,1).value
datas['z']=table.cell(k,2).value
# #重大误区！！！！！
# 原来的代码是这样的saves.append(datas)
# 最后是列表所有元素都一样，因为引用一样
# 所以最后用copy方法添加一个复制，每次循环复制的都不一眼
# 重要
# 重要
# 重要
# 技巧记住！！！！！
# ！！！！
# ！！！！
saves.append(datas.copy())
k+=1
return saves
def write_excel(self):
"""写excel表，写为xls格式"""
new_data=xlwt.Workbook()
sheet1=new_data.add_sheet('结果',cell_overwrite_ok=True)
saves=self.read_excel()
sheet1.write(0,0,'x')
sheet1.write(0, 1, 'y')
sheet1.write(0, 2, 'z')
i=1
for cell in saves:
xd=cell['x']
yd=cell['y']
zd=cell['z']
self.cd.transform(xd,yd,zd)
Xr=self.cd.Xb
print(Xr)
sheet1.write(i,0,Xr[0,0])
sheet1.write(i,1,Xr[1,0])
sheet1.write(i,2,Xr[2,0])
i+=1
file_name=self.file_name_write.get()
#存储在指定路径
print(file_name)
# new_data.save('E:/workspace/BJ54transform_project/Data/test.xls')
new_data.save(file_name)#存储在指定路径
注释比较详细，各位直接看注释就好。
各位小伙伴，永恒的宅邸拉莱耶中，克苏鲁候汝入梦！
扫一扫
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
