#由經(jīng)緯度反算成高斯投影坐標(biāo)（北京54坐標(biāo)系）
#作者：小熊（初學(xué)Python的菜鳥\（￣▽￣）／）

import?math

def?GTXY1（）:

????def?GaussToXY（longitudelatitude）:

????????iPI?=?3.1415926535898/180.0
????????ZoneWide?=?6;?#6度帶寬

????????#54年北京坐標(biāo)系參數(shù)

????????a=6378245.0

????????f=1.0/298.3

????????#80年西安坐標(biāo)系參數(shù)
????????#a=6378140.0;
????????#f=1/298.257;


????????ProjNo?=?（int）（longitude?/?ZoneWide）

????????longitude0?=?ProjNo?*?ZoneWide?+?ZoneWide?/?2?#中央子午線計(jì)算
????????longitude0?=?longitude0?*?iPI?;
????????latitude0=0;

????????longitude1?=?longitude?*?iPI?#經(jīng)度轉(zhuǎn)換為弧度

????????latitude1?=?latitude?*?iPI?#緯度轉(zhuǎn)換為弧度

????????e2=2*f-f*f?#第一偏心率

????????ee=e2*（1.0-e2）

????????NN=a/math.sqrt（1.0-e2*math.sin（latitude1）*math.sin（latitude1））

????????T=math.tan（latitude1）*math.tan（latitude1）

????????C=ee*math.cos（latitude1）*math.cos（latitude1）

????????A=（longitude1-longitude0）*math.cos（latitude1）

????????M=a*（（1-e2/4-3*e2*e2/64-5*e2*e2*e2/256）*latitude1-（3*e2/8+3*e2*e2/32+45*e2*e2*e2/1024）*math.sin（2*latitude1）+（15*e2*e2/256+45*e2*e2*e2/1024）*math.sin（4*latitude1）-（35*e2*e2*e2/3072）*math.sin（6*latitude1））

????????xval?=?NN*（A+（1-T+C）*A*A*A/6+（5-18*T+T*T+14*C-58*ee）*A*A*A*A*A/120）

????????yval?=?M+NN*math.tan（latitude1）*（A*A/2+（5-T+9*C+4*C*C）*A*A*A*A/24+（61-58*T+T*T+270*C-330*ee）*A*A*A*A*A*A/720）

????????X0?=?1000000L*（ProjNo+1）+500000L?

????????Y0?=?0

????????xval?=?xval+X0

????????yval?=?yval+Y0

????????return?xval?yval

????longitude?=?input（‘輸入經(jīng)度:‘）

????print?‘‘

????latitude?=?input（‘輸入緯度:‘）

????print?‘‘

????print?GaussToXY（longitudelatitude）

def?GTXY2（）:

????def?XYToGauss（XY）:

????????iPI?=?3.1415926535898/180.0

????????#54年北京坐標(biāo)系參數(shù)
????????a?=?6378245.0
????????f?=?1.0/298.3

????????#80年西安坐標(biāo)系參數(shù)?
????????#a=6378140.0;
????????#f=1/298.257;

????????ZoneWide?=?6

????????ProjNo?=?（int）（X/1000000L）?#查找?guī)?hào)

????????longitude0?=?（ProjNo-1）?*?ZoneWide?+?ZoneWide?/?2;

????????longitude0?=?longitude0?*?iPI?#中央經(jīng)線

????????X0?=?ProjNo*1000000L+500000L

????????Y0?=?0

????????xval?=?X-X0

????????yval?=?Y-Y0

????????e2?=?2*f-f*f

????????e1?=?（1.0-math.sqrt（1-e2））/（1.0+math.sqrt（1-e2））

????????ee?=?e2/（1-e2）

????????M?=?yval

????????u?=?M/（a*（1-e2/4-3*e2*e2/64-5*e2*e2*e2/256））

????????fai?=?u+（3*e1/2-27*e1*e1*e1/32）*math.sin（2*u）+（21*e1*e1/16-55*e1*e1*e1*e1/32）*math.sin（4*u）+（151*e1*e1*e1/96）*math.sin（6*u）+（1097*e1*e1*e1*e1/512）*math.sin（8*u）

????????C?=?ee*math.cos（fai）*math.cos（fai）

????????T?=?math.tan（fai）*math.tan（fai）

????????NN?=?a/math.sqrt（1.0-e2*math.sin（fai）*math.sin（fai））

????????R?=?a*（1-e2）/math.sqrt（（1-e2*math.sin（fai）*math.sin（fai））*（1-e2*math.sin（fai）*math.sin（fai））*（1-e2*math.sin（fai）*math.sin（fai）））

????????D?=?xval/NN

????????#計(jì)算經(jīng)度（Longitude）?緯度（Latitude）

????????longitude1?=?longitude0+（D-（1+2*T+C）*D*D*D/6+（5-2*C+28*T-3*C*C+8*ee+24*T*T）*D*D*D*D*D/120）/math.cos（fai）;

????????latitude1?=?fai?-（NN*math.tan（fai）/R）*（D*D/2-（5+3*T+

?屬性????????????大小?????日期????時(shí)間???名稱
-----------?---------??----------?-----??----
?????文件????????3614??2017-12-06?16:35??高斯正反算.py