資源簡(jiǎn)介
自己編寫(xiě)的高斯投影正反算Python源碼,經(jīng)驗(yàn)證,精度在0.001以?xún)?nèi),可以實(shí)用,包含北京54和西安80坐標(biāo)系
代碼片段和文件信息
#?高斯坐標(biāo)正反算算法
import?math
P=180/math.pi*3600??#1弧度的秒數(shù)常數(shù)值,1弧度=180/pi=57.29578*3600=206264.80624709636
def?GetGSZS54(BL):?#?北京54克拉索夫斯基橢球參數(shù)??BL?為十進(jìn)制度,非度分秒格式,python中math.sin?和.cos函數(shù)僅支持弧度,需要用math.radians將度數(shù)轉(zhuǎn)換為弧度。
????#北京54橢球高斯正算
????l=(L-Get3DL0(L)[0])*3600/P
????COS2B=math.pow(math.cos(math.radians(B))2)
????N=6399698.902-(21562.267-(108.973-0.612*COS2B)*COS2B)*COS2B
????A0=32140.404-(135.3302-(0.7092-0.0040*COS2B)*COS2B)*COS2B
????A4=(0.25+0.00252*COS2B)*COS2B-0.04166
????A6=(0.166*COS2B-0.084)*COS2B
????A3=(0.3333333+0.001123*COS2B)*COS2B-0.1666667
????A5=0.0083-(0.1667-(0.1968+0.0040*COS2B)*COS2B)*COS2B
????x=6367558.4969*B*3600/P-(A0-(0.5+(A4+A6*math.pow(l2))*math.pow(l2))*math.pow(l2)*N)*math.sin(math.radians(B))*math.cos(math.radians(B))
????y=(1+(A3+A5*math.pow(l2))*math.pow(l2))*l*N*math.cos(math.radians(B))+Get3DL0(L)[1]*1000000+500000??#坐標(biāo)加3度帶帶號(hào),y坐標(biāo)西移500KM
????return(xy)
def?GetGSZS80(BL):?????#?西安80??國(guó)際1975國(guó)際橢球參數(shù)
????#西安80國(guó)際1975橢球高斯正算
????l=(L-Get3DL0(L)[0])*3600/P
????COS2B=math.pow(math.cos(math.radians(B))2)
????N=6399596.652-(21565.045-(108.996-0.603*COS2B)*COS2B)*COS2B
????A0=32144.5189-(135.3646-(0.7034-0.0041*COS2B)*COS2B)*COS2B
????A4=(0.25+0.00253*COS2B)*COS2B-0.04167
????A6=(0.167*COS2B-0.083)*COS2B
????A3=(0.3333333+0.001123*COS2B)*COS2B-0.1666667
????A5=0.00878-(0.1702-0.20382*COS2B)*COS2B
????x=6367452.1328*B*3600/P-(A0-(0.5+(A4+A6*math.pow(l2))*math.pow(l2))*math.pow(l2)*N)*math.cos(math.radians(B))*math.sin(math.radians(B))
????y=(1+(A3+A5*math.pow(l2))*math.pow(l2))*l*N*math.cos(math.radians(B))+Get3DL0(L)[1]*1000000+500000??#坐標(biāo)加3度帶帶號(hào),y坐標(biāo)西移500KM
????return(xy)
def?GetGSFS54(xyL0):
?????#北京54橢球高斯反算
????y=y-L0/3*1000000-500000??#去掉帶號(hào)和y坐標(biāo)恢復(fù)東移500KM
????BT=x/6367558.4969??#??BT=x/6367558.4969*P?得到是秒數(shù),不乘P直接得到弧度數(shù)值
????bt=x/6367558.496
評(píng)論
共有 條評(píng)論