高斯投影算法

① 高斯平面直角坐标的投影方法

高斯投影的方法是将地球按经线划分为带，称为投影带。投影是从首子午线开始的，分6°带和3°两种。每隔6°划分一带的叫6°带，每隔3°划分一带的叫3°带。我国领土位于东经72°∽136°之间，共包括了11个6°带，即13∽23带；22个3°投影带即24∽45带。
设想一个平面卷成横圆柱套在地球外，如图1－5（a)所示 。通过高斯投影，将中央子午线的投影作为
纵坐标轴，用x表示，将赤道的投影作横坐标轴，用y表示，两轴的交点作为坐标原点，由此构成的平面直角坐标系称为高斯平面直角坐标系，如图1－5（b) 所示。每一个投影带都有一个独立的高斯平面直角坐标系，区分各带坐标系则利用相应投影带的带号。在每一个投影带内，y坐标值都有正有负，这对于计算和使用都不方便，为了使y坐标都为正值，故将纵坐标轴向西平移500㎞，并在y坐标前加上投影带的带号。 6°带投影是从英国格林尼治子午线开始，自西向东，每隔经差6°分为一带，将地球分为60个带，其编号分别为1，2，3，…60。任意带的中央子午线经度为Lo，它与投影带号N的关系如下所示：
Lo=(6N-3°)
式中：N———6°带的带号
离中央子午线越远，长度变形越大，在要求较小的投影变形时，可采用3°投影带。3°带是在6°带的基础上划分的，如图所示。每3°为一带，从东经1°30′开始，共120带，其中央子午线在奇数带时与6°带的中央子午线重合，每带的中央子午线可用下面的工式计算：
Lo=3N′
式中：N′——3°带的带号。
为了避免y坐标出现负值，3°带的坐标原点同6°带一样，向西移动500㎞，并在y坐标前加3°带的带号。

② 高斯-克吕格投影理解

你想得不错，确实有你说的这种投影的应用，那就是大名鼎鼎的“墨卡托投影”。其实高斯—克吕格投影以及UTM投影都可算作是墨卡托投影的变种。高斯投影是“等角横切圆柱投影”，墨卡托投影时"等角正切圆柱投影”，UTM与墨卡托更类似，只不过用于军事（美国）。高斯投影没有角度变形——这是最关键的，在长度和面积上变形也很小，中央经线肯定无变形，自中央经线向东向西，变形逐渐增加。墨卡托投影的长度和面积变形明显，但标准纬线无变形，从标准纬线向两极变形逐渐增大，但因为它具有各个方向均等扩大的特性，保持了方向和相互位置关系的正确——这是最重要的。因此常用作航海图和航空图——也就是说绝对方位他不准，但是相对位置他最准。

③ 用通俗的方法讲解高斯投影和中央子午线

中央子午线" 英文对照
central meridian;

"中央子午线" 在工具书中的解释

1、又称“中央经线”。位于投影带中央的子午线。中央经线一般为直线,其他经线分布在它的两侧呈弧线。高斯投影带的中央子午线即为一条直线,其长度不变。在六度带中,它的经度为L=6°×n—3°,n为六度带的带号。(参看高斯投影分带)。在其他小比例尺地图投影中,中央经线也为直线,多通过所表示的主要地区。

"中央子午线" 在学术文献中的解释

1、每一个投影带的中间一条子午线称为中央子午线,其经度为6°n-3°(n为中央子午线的编号).而3°带的中央子午线经度,则是按经3°为一带,其经度是3°n(n为3°带的编号)
文献来源

2、高斯投影是设想将截面为椭圆的一个圆柱面横套在旋转椭球外面（图A）并与旋转椭球面上某一条子午线（NOS）相切同时使圆柱的轴位于赤道面内并通过椭球中心相切的子午线称为中央子午线.然后将中央子午线附近的旋转椭球面上的点、线投影到横圆柱面上如将旋转椭球面上的M点投影到横圆柱面上得m点再顺着中央子午线将圆柱面剪开展成平面如图B所示这个平面称为高斯投影平面

===============================================================
1 前言
在地球物理勘察中，常常要提供测区、有利成矿区的面积，常用到面积量算。在小范围内的测量中，通常我们都是先测出区域边界点的平面坐标，再按封闭的区域计算面积，而在大面积的勘察中，这样计算面积就不合适了，因为，测量上是以参考椭球面作为地球的参考面的，表现在该面上的地面图形是曲面上的图形，由于地球半径很大，当测区面积较小时，可以把参考椭球面上的曲面当水平面看待，采用水平投影的方法即可，但测区面积较大时，曲面就不能当成平面看待了。
2 高斯投影的概念
高斯投影又称横轴椭圆柱等角投影，是德国测量学家高斯于1825～1830年首先提出的。实际上，直到1912年，由德国另一位测量学家克吕格推导出实用的坐标投影公式后，这种投影才得到推广，所以该投影又称高斯-克吕格投影。高斯投影就是设想一个横椭园柱面作为投影面的分带正形投影。如图所示，使椭圆柱的轴通过旋转椭圆体中心且与旋转椭圆体的长轴重合(即与赤道面重合)，同时使椭圆柱面与旋转椭圆体投影带的中央子午线相切，在保持等角条件下，用数学的方法，将旋转椭圆体投影带上的点、线投影到横椭圆柱面上。如旋转椭圆体面上的A点投影到椭圆柱面上为a点，赤道面与椭圆柱面的交线EE为赤道的投影。 投影后，依过极点的母线（平行于椭圆柱轴的椭圆柱面上的直线）将椭圆柱面切开，并展成平面M，如图所示，该平面叫高斯投影平面。
高斯投影有以下三条规律：
（1）旋转椭球体面上两极间经差相等的投影带的子午线，除中央子午线投影后长度不变且为直线外，其余子午线投影后为凹向中央子午线的曲线，并以中央子午线为对称轴，其长度大于投影前长度，且离中央子午线愈远变形愈大，因此投影带上的线段除位于中央子午线上者外，其投影后长度均较实地长度增长。据推算，离中央子午线300公里时投影长度变形为1/900，而且由一点出发各方向的变形是成比例的。
（2）投影后的纬圈除赤道为直线外，其余均为凸向赤道的曲线，并以赤道为对称轴。
（3）经线与纬圈投影以后仍然保持正交。中央子午线与赤道投影后，变成互相垂直的直线。
3 高斯投影对面积测量的影响
高斯投影对面积测量的影响主要是通过对长度的影响引起的，我们进行地面测量实际上就是对大地椭球进行测量，根据高斯投影的概念我们知道，将椭球面投影到高斯平面上会产生一定的误差，误差的形成主要有两个方面的原因：
（1）远离中央子午线引起的误差
高斯投影对测区面积测量的影响与测区到中央子午线的距离有关，离中央子午线越远误差就越大，我们以纬度为35°，离中央子午线不同的距离，经差纬差都是1′的区块为例，其在椭球面上的面积是2813026.48m2，在高斯平面上因远离中央子午线有不同的面积，如下表
纬差 10′ 30′ 1° 2° 3°
面积(m2) 2813044.16 2813175.69 2813613.67 2815356.60 2818257.47
面积比 1.000006285 1.000053 1.000209 1.000828 1.00186
表中面积比为高斯平面上面积与椭球面上面积的比，由表中可以看出，离中央子午线越远面积增大得越多。
（2）投影面高程引起的误差
将椭球面投影到不同高度的椭圆柱面上引起的误差也是不同的，投影面高程越高误差也越大，我们选择椭球面上四个点，经纬度分别为 （111°00′30〃， 35°00′00〃），（110°59′30〃， 35°00′00〃），（110°59′30〃， 35°01′00〃），（111°00′30〃 ，35°01′00〃），通过计算我们得出其在椭球面上的面积为2813026.48m2
下表列举了不同投影面高程对面积的影响
投影面高程 100m 200m 500m 1000m 2000m 5000m
面积(m2) 2813114.62 2813202.76 2813467.2 2813907.95 2814789.57 2817435.25
面积比 1.00003133 1.00006267 1.00015667 1.00031335 1.00062676 1.00156727
表中面积比为高斯面上面积与椭球面上面积的比。通过研究我们得出结论：远离中央子午线、投影面高程使面积增大。
注：以上各表均采用武汉中地信息工程有限公司研制的mapgis地理信息系统进行计算。
4 正确计算测区面积的方法
（1）如果采用高斯投影的方法来计算面积，要根据本区的平均海拔高程选择合适的投影面高程，选择三度带或更小的分带来进行投影，若遇投影变形较大影响到测量精度时，还可以采用独立坐标系统的任意带投影，这样就减小了投影面高程和中央子午线对面积的影响。
（2）采用其它投影方式如亚尔勃斯等积圆锥投影，这种投影是按等面积条件，将地球上的经纬线投影到剖于地球某两条平行圈的圆锥面上，沿一条母线将圆锥面展开成平面。选用这种投影可以保证投影面面积与椭球面面积相等。当然通常只在计算面积时选用这种投影，其它测量工作还可以用高斯投影来进行。
（3）我们在探矿权测量中引入了区块的概念，即将矿区分成若干很小的区块，比如经差纬差都是15〃的小区块，因为在相同纬度上，相同每一个小区块的面积是相等的，这样我们只要直接计算出一条经线上每隔15〃一个小区块在地球椭球面上的面积，再计算出不同纬度区块的数量就可以计算出矿区的面积了。

④ 哪里能够找到高斯投影反算公式

这是一个高斯反算的算法，你可以参考以下，这个程序是针对克拉索夫斯基椭球的，转成其他的椭球体的坐标，只需要换一下参数就可以了。
#define PI 3.1415926535
main()
{
float P;
float x,y,z;
float a,b,bf,Bf,B,l1,l2,L,L0;
float b2,b3,b4,b5;
float N,i,j,k;
int B1,B2,L1,L2;float B3,L3;
P=(180.0/PI)*3600.0;
printf("please enter x,y:\n");
scanf("%f,%f",&x,&y);
y-=500000;
a=x*P/6367558.4969;
b=(a*PI/180.0)/3600.0;
i=cos(b)*cos(b);
bf=a+(50221746+(293622+(2350+22*i)*i)*i)*sin(b)*cos(b)*P*(1E-10);
Bf=(bf*PI/180.0)/3600.0;
j=cos(Bf)*cos(Bf);
N=6399698.902-(21562.267-(108.973-0.612*j)*j)*j;
b2=(0.5+0.003369*j)*sin(Bf)*cos(Bf);
b3=0.333333-(0.166667-0.001123*j)*j;
b4=0.25+(0.16161+0.00562*j)*j;
b5=0.2-(0.1667-0.0088*j)*j;
z=(y/N)/cos(Bf);
k=z*z;
printf("Please enter L0:\n");
scanf("%f",&L0);
B=bf-(1-(b4-0.12*k)*k)*k*b2*P;
B=B/3600.0;
B1=(int)B;
B2=(int)((B-B1)*60);
B3=((B-B1)*60-(int)((B-B1)*60))*60.;
l1=(1-(b3-b5*k)*k)*z*P;
L=(l1/3600)+L0;
L1=(int)L;
L2=(int)((L-L1)*60);
L3=((L-L1)*60-(int)((L-L1)*60))*60.;
}
}

⑤ 计算校正法

(1)收集测区GPS控制点

一般情况使用的基础地形图是北京54坐标系或西安80坐标系，这样就应收集该区域已知点的北京54坐标或西安80坐标，可以从当地测绘部门收集到，一般收集的坐标是北京54坐标系或西安80坐标系的平面坐标(X，Y，H)，这里的H通常是重力高程即海拔高度，所以需要进行换算获取椭球高程(如果没有换算参数也可以直接使用，一般不会影响平面坐标)。首先将获取的平面坐标系坐标(X，Y)反算求出相应大地坐标(B，L)，这里的X、Y是高斯投影坐标，计算使用高斯投影的反变换公式，这里不再列举，也可以使用软件进行计算如MapGIS、ERDAS都具有这项功能，然后由大地坐标计算空间直角坐标系坐标。

通常情况下至少要收集三个控制点，且这三个控制点不应在一条直线上。然后收集这些控制点的WGS84的大地坐标(B，L，H)或直接测量出这些控制点的WGS84的大地坐标(B，L，H)，并换算成WGS84空间直角坐标系坐标(X，Y，Z)，这样就可以计算dX、dY、dZ。

(2)计算坐标转换的其他三个参数

dX、dY、dZ三个参数就是根据前面求出的控制点的WGS84与北京54或西安80大地坐标相对应的空间直角坐标系的X、Y、Z坐标的差值，当测区具有多个控制点时，取相应平均值作为参数值dX、dY、dZ。

数字地质填图技术实习教程

(3)计算坐标转换参数的具体方法

由以上论述可以知道计算三个参数的理论方法，那么如何进行实际的计算呢?由高斯投影平面坐标向大地坐标的转换，GIS软件或者遥感软件都支持比如MapGIS、ERDAS，都可以完成相关的计算，所以这一步可以采用现有的软件完成。

由大地坐标向地心直角坐标的转换，现有的GIS或者遥感软件鲜有支持者，所以需要进行转换。为了方便起见，将转换的公式再次列于下面。

设椭球长半轴为a，扁率为f，H为相对椭球面的高度，那么大地坐标(B，L)处地心直角坐标系的计算公式为

数字地质填图技术实习教程

式中:N为纬度B处的卯酉圈曲率半径，N=a/(1－e²sin²B)^1/2;e为椭球第一偏心率，e²=2f－f²。

大地坐标向地心直角坐标的转换有两种方法:

第一种方法是使用Excel表格进行计算，将上述公式输入Excel中，输入数据就可以实现计算的自动化。

第二种方法是自己编写计算程序进行计算，可以使用任何计算机语言进行编写，下面列出了使用C++语言编写的计算转换参数程序的源代码，需要指出的是，这一程序可以有许多不同的编写方法，这里只是列出了其中一种写法以供参考，并且在这个程序中并没有使用C++的标准库，这主要是为了简化程序，使用标准类库会带来更高的效率和稳定性，同时此程序没有进行任何异常或者错误的处理，请读者自己完成，以保证程序的正确性。

∥自定义结构体:ReferenceEllipseCoor，别名RECoor

∥字段:B、L、H

∥目的:存储大地坐标数据，作为指针参数传递

∥大地坐标系结构

{

doubleB;∥大地纬度

doubleL;∥大地经度

doubleH;∥椭球高

}RECoor;

∥函数名:

∥目的:计算坐标系转换用da、df、dX、dY、dZ五参数

∥参数As:原椭球体长半轴

∥参数At:目标椭球体长半轴

∥参数Fs:原椭球体扁率

∥参数Ft:目标椭球体扁率

∥参数pSPoints:原椭球体上的大地坐标，RECoor指针类型

∥参数pTPoints:目标椭球体上的大地坐标，RECoor指针类型

∥参数n:转换用坐标点对个数

∥参数da:引用类型返回转换参数之一DA

∥参数df:引用类型返回转换参数之一DF

∥参数dX:引用类型返回转换参数之一DX

∥参数 dY: 引用类型返回转换参数之一 DY

∥参数 dZ: 引用类型返回转换参数之一 DZ

void (double As，double At，

double Fs，double Ft，RECoor * pSPoints，RECoor * pTPoints，

int n，double＆ dA，double＆ dF，double＆ dX，double＆ dY，double＆ dZ)

{

double dXSum = 0;

double dYSum = 0;

double dZSum = 0;

double Xs;

double Ys;

double Zs;

double Xt;

double Yt;double Zt;

double esq = 2* Fs － Fs* Fs;

double etq = 2* Ft － Ft* Ft;

for(int i = 0; i < n; i + + )

{

double B = pSPoints ［i］ . B;

double L = pSPoints ［i］ . L;

double H = pSPoints ［i］ . H;

double SinB = : : sin(B);

double CosB = : : cos(B);

double SinB2 = SinB* SinB;

double SinL = : : sin(L);

double CosL = : : cos(L);

double N = As/: : sqrt(1 － esq* SinB2);Xs =(N + H)* CosB* CosL;

Ys =(N + H)* CosB* SinL;

Zs =(N*(1 － esq)+ H)* SinB;

B = pTPoints ［i］ . B;

L = pTPoints ［i］ . L;

H = pTPoints ［i］ . H;

SinB = : : sin(B);

CosB = : : cos(B);

SinB2 = SinB* SinB;

SinL = : : sin(L);

CosL = : : cos(L);

N = At/: : sqrt(1 － etq* SinB2);

Xt =(N + H)* CosB* CosL;

Yt =(N + H)* CosB* SinL;

Zt =(N*(1 － etq)+ H)* SinB;

dXSum + =(Xs － Xt);

dYSum + =(Ys － Yt);

dZSum + =(Zs － Zt);

}

dA = As-At;

dF = Fs-Ft;

dX = dXSum/n;

dY = dYSum/n;

dZ = dZSum/n;

}

⑥ 求高斯投影中的中央经线

1∶2.5万及1∶5万的地形图采用6度分带投影。

中央子午线为102度，带号为17，(6度带投影) 102/6=17。

======
原理如下：
1∶2.5万及1∶5万的地形图采用6度分带投影，即经差为6度，从零度子午线开始，自西向东每个经差6度为一投影带，全球共分60个带，用1，2，3，4，5，……表示．即东经0～6度为第一带，其中央经线的经度为东经3度，东经6～12度为第二带，其中央经线的经度为9度……。

1∶1万的地形图采用3度分带，从东经1.5度的经线开始，每隔3度为一带，用1，2，3，……表示，全球共划分120个投影带，即东经1.5～ 4.5度为第1带，其中央经线的经度为东经3度，东经4.5～7.5度为第2带，其中央经线的经度为东经6度．

如何计算当地的中央子午线？
当地中央子午线决定于当地的直角坐标系统，首先确定您的直角坐标系统是3度带还是6度带投影
---------------------------------------------------------------------------------------------
公式推算:
6度带中央子午线计算公式：当地经度/6=N；中央子午线L=6 * N （带号）
当没有除尽，N有余数时， 中央子午线L=6*N - 3
---------------------------------------------------------------------------------------------
3度带中央子午线计算公式： 当地经度/3=N；中央子午线L=3 X N
---------------------------------------------------------------------------------------------
我国的经度范围西起 73°东至135°，可分成
六度带十一个（13号带—23号带），各带中央经线依次为（75°、81°、……123°、129°、135°）；

三度带二十二个（24号带—45号带）。各带中央经线依次为（72°、75°、……132°、135°）；

六度带可用于中小比例尺（如 1：250000）测图，三度带可用于大比例尺（如 1：10000）测图，城建坐标多采用三度带的高斯投影
---------------------------------------------------------------------------------------------
如何判断投影坐标是3度带坐标还是6度带坐标
如(4231898,21655933)其中21即为带号，同样所定义的东伪偏移值也需要加上带号，如21带的东伪偏移值为21500000米。假如你的工作区经度在120度至126度范围，则该坐标系为6度带坐标系，该带的中央经度为123度。
如(2949320,36353822)其中36即为带号，已知该地点位于贵阳市附近，而从地图上我们看到贵阳大概的经度是东经108度左右，因此可以36*3=108，所以该坐标系为3度带坐标系，该带的中央经度为108度。而不可能为6度带：36*6=216。高斯投影是正形投影的一种，同一坐标系中的高斯投影换带计算公式是根据正形投影原理推导出的两个高斯坐标系间的显函数式。在同一大地坐标系中(例如1954北京坐标系或1980西安坐标系)，如果两个高斯坐标系只是主子午线的经度不同，那么显函数式前的系数可以根据坐标系使用的椭球元素和主子午线经度唯一确定。但如果两个高斯坐标系除了主子午线的经度不同以外，还存在其他线性系，则将线性变换公式代入换带计算的显函数式中，仍然可以得到严密的坐标变换公式。此时显函数式前的系数等价于使用两个坐标系主子午线的经度和线性变换参数联合求解得到的，可以唯一确定。

⑦ 高斯投影带3度分带法为什么从东经1°30′起开始

应该是为保证3°分带和6°分带的中央经线一致

⑧ 高斯投影是通过什么方法将地球椭球面上的图形展绘到平面上的

地球表面的距离，经过高程归化改正，归化到地拖椭球面上，然后通过高斯投影改正，计算到高斯平面上。

⑨ 2．设我国某地一点的高斯坐标为： ＝3234567.89m， ＝38432109.87m。问该点是属于几度带高斯投影坐标

你这个坐标是3度带的。
3度带的中央子午线是114,算法是Y的前2位乘以3,你这个坐标是3度带的.我国的经度范围是西起73,东至135，6度带范围为13-23,一共11个，3度带范围为24-45，共22个。

⑩ 年北京坐标系与西安坐标系的转换方法

在矿业权实地核查准备工作阶段，收集到的地质、测绘等相关资料、图件和矿业权登记数据中，所涉及的地理数据可能是不同大地坐标系下的坐标数据。从实际情况来看，矿业权拐点坐标大多采用的是1954年北京坐标系，矿区已有的测量控制点和测量资料多数采用的也是1954年北京坐标系。本次矿业权实地核查测量工作采用的是1980西安坐标系，在实地测量和数据整理中涉及1954年北京坐标系与1980西安坐标系的转换。下面简要介绍二者之间转换的理论与方法。

（一）高斯投影正算和反算

将大地坐标换算为平面直角坐标，叫做高斯投影正算，是在同一椭球中进行，不存在误差。其常用量定义和公式如下：

a为椭球长半轴

b为椭球短半轴

f为椭球扁率

e为第一偏心率

e＇为第二偏心率

全国矿业权实地核查技术方法指南研究

B为纬度，单位为弧度

全国矿业权实地核查技术方法指南研究

M为子午圈曲率半径

N为卯酉圈曲率半径

子午线弧长X

设有子午线上两点p1和p2,p1在赤道上，p2的纬度为B,p1、p2间的子午线弧长X计算公式：

全国矿业权实地核查技术方法指南研究

例如，1980西安坐标系a=6378140,e²=0.006694385,A'＝1.005052506,B'＝0.002531556209,C'＝2.656901555E-06,D'＝3.470075599E-09,E'＝4.916542167E-12,F '＝7.263137253E-15,G'＝1.074009912E-17以B=30°弧度值0.5235987756为例，在Y=0时算得X=3320114.946。

当Y≠0,l≠0时则需要采用下列积分和逐次趋近的方法。

（1）高斯正算公式（利用点的经纬度计算XY坐标）

全国矿业权实地核查技术方法指南研究

（2）高斯反算公式（利用点的XY坐标计算经纬度）

全国矿业权实地核查技术方法指南研究

（3）底点纬度B_f迭代公式

全国矿业权实地核查技术方法指南研究

直到B_i-1-B_i小于某一个指定数值，即可停止迭代。

式中

全国矿业权实地核查技术方法指南研究

国家测绘局经过改进，将7个系数改为5个算出各椭球的值，采用公式如下：

（1）高斯投影正算（B,L→x,y）

全国矿业权实地核查技术方法指南研究

式中：X₀=C₀B-cosB（C₁sinB+C₂sin²B+C₂sin⁵B+C₄sin⁷B）

m₀=lcosB

l=L－中央子午线经度值（弧度）

L,B为该点的经纬度值。

全国矿业权实地核查技术方法指南研究

式中：t=tanB，η²=e'²cos²B,

C,C₀,C₁,C₂,C₃,C₄,e²为椭球常数

（2）高斯投影反算（x,y→B,L）

全国矿业权实地核查技术方法指南研究

式中：t=tanB_f，η²=e'²cos²B_f，

，K₁K₂,K₃,K₄为椭球常数。

各坐标系椭球常数如表4-1。

表4-1 各大地坐标系椭球常数

国家测绘局采用的公式编程更加容易，高斯投影的正算、反算因为是在同一椭球下进行，公式是严密的，不存在误差，电算操作非常方便。现在网上很多软件有这种功能。度、分、秒输入使用小数形式，小数点前面是度，小数点后前两位为分，后两位为秒，再后面为秒的十进制小数。如25.23451124其值为25°23′45.1124″，正反算已经成了非常简单的事。高斯正算、反算必须考虑到椭球参数，椭球不同结果是不同的。必须考虑到中央子午线位置。因为各带中都有重复点，本次实地核查要求使用3度带，所有Y坐标必须带有3°带的带号，不允许使用独立坐标系或假定坐标系。

（二）参心坐标与空间直角坐标的关系

空间直角坐标X、Y、Z与大地坐标B、L、H间的关系表示如下：

全国矿业权实地核查技术方法指南研究

大地坐标B、L、H 与空间直角坐标X、Y、Z间的关系表示如下：

全国矿业权实地核查技术方法指南研究

式中

在转换中对于不知道椭球高的控制点可将控制点的大地高置为0，放在椭球面上计算，三维就变成二维，其效果更好。

（三）坐标系统转换

1954年北京坐标系与1980西安坐标系的转换通常有两种方法：四参数转换法和七参数转换法。

1.四参数转换法

所谓四参数转换是两个平移参数，一个旋转参数，一个尺度比。不考虑什么椭球，在小范围内按平面坐标直接平移、旋转、缩放。最少条件是两个公共点，多公共点时可以使用最小二乘法，删除残差大的点。这在区域面积小的情况下是可以的，一般不宜超过40平方千米。四参数转换模型如下：

x₂＝Δx+x₁（1+m）cosa-y（1+m）sina

y₂＝Δx+x₁（1+m）sina-y（1+m）cosa

2.七参数转换法

该方法适用于椭球间的坐标转换。其实质是原椭球空间直角坐标（X₁,Y₁,Z₁）与新椭球空间直角坐标（X₂,Y₂,Z₂）间的转换。椭球间的坐标转换至少需要3个公共点，解算七参数。转换公式采用的是布尔莎公式，法方程的解算采用高斯消元法。高斯消元法，是线性代数中的一个算法，可用来为线性方程组求解，求出矩阵的秩，以及求出可逆方阵的逆矩阵。当用于一个矩阵时，高斯消元法会产生出一个“行梯阵式”。高斯消元法可以用在电脑中来解决数千条等式及未知数。迭代法较消元法的残差大。

椭球间的坐标转换适用基于椭球的参心（地心）坐标系间的转换，而不适用于基于平面的独立坐标系间以及独立坐标系和参心（地心）坐标系间的转换。基于椭球的坐标转换中（七参数），椭球→椭球的转换实际上是在空间直角坐标系中完成的。完整的变换过程如下（以“平面→平面”为例）：（x₁,y₁,H₁）→（B₁,L₁,H₁）→（X₁,Y₁,Z₁）→（X₂,Y₂,Z₂）→（B₂,L₂,H₂）→（x₂,y₂,H₂）。首先把直角坐标系下的直角坐标，原公共点中的1954年北京坐标转换成2000国家大地经纬度坐标，再转换为1954年北京坐标系的参心坐标，公共点的1980西安坐标做同样转换。利用两个椭球的参心（地心）坐标求得转换参数，利用该参数直接将1954年北京坐标系下的坐标转换成1980西安坐标系下的坐标。在上述过程中，高程H₁、H₂是大地高（椭球高）。大地高＝正常高＋测区高程异常。如果不需要转换高程的话，可以将高程和高程异常全部置为0。不可将1954年北京坐标系坐标所带的正常高直接代入。

七参数的转换模型如下：

（1）七参数转换模型

全国矿业权实地核查技术方法指南研究

式中：ΔB，ΔL为同一点位在两个坐标系下的纬度差、经度差（弧度）；

a，Δf为椭球长半轴差（米）、扁率差（无量纲）；

X，ΔY，ΔZ为平移参数（米）；

ε_x，ε_y，ε_z为旋转参数（弧度）；

m为尺度参数（无量纲）。

最少3个公共点可以解求出七个参数。

（2）三维七参数转换模型

全国矿业权实地核查技术方法指南研究

全国矿业权实地核查技术方法指南研究

式中：ΔB，ΔL，ΔH为同一点位在两个坐标系下的纬度差（弧度）、经度差（弧度）、大地高差（米）；

ρ为一个弧度的秒值，180×3600/π弧度/秒；

a为椭球长半轴差（米）；

f为扁率差（无量纲）；

X，ΔY，ΔZ为平移参数（米）；

εx，ε_y，ε_z为旋转参数（弧度）；

m为尺度参数（无量纲）。

最少3个公共点可以解求出七个参数。

七参数适用于整个测区的转换，面积小于2000平方千米的可以一次转换完成，面积大的可以分区转换，各分区之间应选公共点，以保证数据的接边精度。关于残差，国家规定以1∶2000图为例，残差为图上0.1毫米即实地20厘米，超过3倍中误差的点删除。为了保证矿业权矿界拐点转化的精度，本次矿业权实地核查规定残差超过实地0.1米一般不宜使用，实际上比国家规定的精度严，相当于国家规定的1/6。

（四）利用坐标转换软件进行坐标转换

以上介绍了1954年北京坐标系和1980西安坐标系转换的理论，在实际转换时可以采用相关的软件来完成。目前，市场上有多种坐标转换软件可供选择。在选择软件时，应注意部分软件转换的精度可能达不到本次矿业权实地核查的要求。下面以经天测绘技术公司开发的测量计算工具包软件V4.05为例，介绍坐标转换方法。

该软件界面如图4-3。该软件可以进行高斯正算、高斯反算、坐标换带、椭球间的转换，可以批量导入，可以保存数据、保存公共点，包括了坐标转换所需的相关计算功能。另外，该软件还能实现2000国家大地坐标系与1954年北京坐标系、1980西安坐标系、WGS-84坐标系以及独立坐标系的转换。

图4-3 经天测绘技术公司开发的测量计算工具包软件界面

坐标系统变换，可以采用平面坐标转换中的多公共点相似变换和椭球坐标转换。小面积可以采用多公共点相似变换。限制在400平方千米左右，不超过1 幅1∶50000图。它与中央子午线无关、高程需要置为0，计算参数的输入文件为文本文件，格式为：

点号，原X 坐标，原Y坐标，新X 坐标，新Y坐标

需要转换的输入文件格式为：

点号，原X 坐标，原Y坐标

参数计算点数不超过30个，文件可以导入，公共点可以保存，参数也可以保存。转换坐标可以导入，转换后的坐标可以保存。需要注意的是，转换坐标的位数与计算参数的坐标位数应一致。计算参数不使用带号，转换后坐标也没有带号。图4-4中的算例X舍去前4位，Y舍去前3位。

图4-4 多公共点平面相似变换窗口

面积较大的测区应使用7参数转换。在椭球间坐标转换开关下，有平面－平面、大地－平面、平面－大地、大地－大地4个子开关。对于采矿权，可使用平面－平面；对于探矿权，使用大地－大地，小数后位数较多，根据需要可将尾部删去。输入文件的格式与上述相同，需要输入中央子午线，Y坐标不加带号，在不知道1954年北京坐标、1980西安坐标的椭球高的情况下，可在高程栏输入0，测区高程异常输入0，探矿权是大地坐标格式，小数点前3位为°，后2位为′，3、4位为″，后面为十进制的秒的小数，如108°33′15″8563，输入108.33158563，由于控制点坐标是X、Y格式，可用高斯投影反算将控制点变为大地坐标格式。或是使用高斯坐标正算把探矿权登记坐标转换为直角坐标，计算完成后再使用高斯坐标反算将1980西安坐标转换为2000国家坐标。图4-5表示一个县的采矿权转换过程，Y坐标略去了前3位数。

图4-5 椭球间平面坐标转换窗口

需要注意的是，该软件没有采用软件狗加密，但需要注册才能用，采用机器码注册，一个软件只能装一台计算机专用。