聊聊GIS中那些坐标系

转载请声明到标题。

B站/博客园/CSDN/知乎:@岭南灯火

很开心能跃居百度关键词“GIS 坐标系”第一位。

已重写本篇,希望读者们移步:

[聊聊GIS中的坐标系|再版]

从第一次上地图学的课开始,对GIS最基本的地图坐标系统就很迷。也难怪,我那时候并不是GIS专业的学生,仅仅是一门开卷考试的专业选修课,就没怎么在意。

等我真正接触到了各种空间数据产品,我才知道万里长征第一步就是:处理坐标系统。

想必各位从业人员多多少少都会听说过几个名词,可能有那么点印象吧。比如,高斯克吕格,北京54,西安80,WGS84,投影坐标系统等等。

今天就从头说起,讲讲那些坐标系统的事情。

惯例,给个目录:

  1. 经纬度与GCS(Geographic Coordinate System, 地理坐标系统)
  2. 平面坐标与PCS(Projection Coordinate System, 投影坐标系统)
  3. GCS和PCS的转化问题(三参数与七参数问题)
  4. 火星坐标问题
  1. 在第一部分,我介绍一下以经纬度为准的地理坐标系统,也顺带提及一下我国的高程坐标系。主要涉及的内容有:大地水准面问题,椭球问题,常见的GCS(如北京54,西安80,CGCS2000,WGS84等),让大家看到GIS数据中的GCS马上就能知道这是什么东西。
  2. 在第二部分,我介绍一下以平面直角坐标系为量度的投影坐标系统。主要涉及的内容有:PCS与GCS的关系,我国常见的PCS(高斯克吕格、兰伯特/Lambert、阿尔伯斯Albers、墨卡托Mercator、通用横轴墨卡托UTM、网络墨卡托Web Mercator)。
  3. 在第三部分,是实际操作过程中遇到的种种问题,如投影不对会出现什么情况、如何转换GCS、如何切换PCS(重投影问题)等问题,涉及一些数学转换的思维,需要有一定的空间想象能力。
  4. 在第四部分,我简单介绍一下所谓的火星坐标。

那么我们开始吧!

1. 经纬度与GCS

天气预报也好,火箭发射也罢,地震、火山等事故发生时,电视台总会说东经XX度,北纬YY度。这个经纬度中学地理就学过了,我就不细说了。

我从如何描述地球说起。

1.1 凹凸不平的地球

谁都知道地球表面不平坦,它甚至大概形状都不是一个正球体,是一个南北两极稍扁赤道略胖的胖子,胖度大概是20km,在外太空几乎看不出来的,这也可能和星球长期受到潮汐引力、太阳引力以及自身旋转的向心力有关。这里不是地球科学,就不再深究了。

为了能让地球出现在数学家的公式里,我们曾经走过了2个阶段:用平静的海面描述地球——用虚拟的旋转椭球面描述地球表面。

这里也不是地图学,再深入下去其实还有似大地水准面等概念。就挑重点讲。

“假设地球表面都是水,当海平面风平浪静没有波澜起伏时,这个面就是大地水准面。”大家应该知道,在太空失重的环境下,水相对静止状态是个正球体,那么肯定很多人就认为,大地水准面就是个正球面。不是的,还需要考虑一个问题:地球各处的引力不同。引力不同,就会那儿高一些,这儿低一些,尽管这些微小的差距肉眼难以观测出来,可能隔了好几千米才会相差几厘米。所以,在局部可能看起来是个球面,但是整体却不是。显然,用大地水准面来进行数学计算,显然是不合适的,至少在数学家眼中,认为这不可靠。

所以找到一个旋转椭球面就成了地理学家和数学家的问题。(注意区分椭球面和旋转椭球面这两个数学概念,在GCS中都是旋转椭球面)

给出旋转椭球面的标准方程:

(x2+y2)/a2+z2/b2=1

其中x和y的参数相同,均为a,这就代表一个绕z轴旋转的椭圆形成的椭球体。不妨设z轴是地球自转轴,那么这个方程就如下图是一个椭球体,其中赤道是个圆。

这样,有了标准的数学表达式,把数据代入公式计算也就不是什么难事了。

由此我们可以下定义,GIS坐标系中的椭球,如果加上高程系,在其内涵上就是GCS(地理坐标系统)。其度量单位就是度分秒。

描述一个旋转椭球面所需的参数是方程中的a和b,a即赤道半径,b即极半径,f=(a-b)/a称为扁率。

与之对应的还有一个问题:就是坐标中心的问题。(地球的中心在哪里?)

【注】十九世纪发现赤道也是一个椭圆,故地球实际应以普通椭球面表示,但是由于各种原因以及可以忽略的精度内,一直沿用旋转椭球体作为GCS。

1.2 参心坐标系、地心坐标系

上过中学物理的人知道,物体均有其质心,处处密度相等的物体的质心在其几何中心。所以,地球只有一个质心,只是测不测的精确的问题而已。由地球的唯一性和客观存在,以地球质心为旋转椭球面的中心的坐标系,叫地心坐标系,且唯一。当然,由于a、b两个值的不同,就有多种表达方式,例如,CGCS2000系,WGS84系等,这些后面再谈。

【注】地心坐标系又名协议地球坐标系,与GPS中的瞬时地球坐标系要对应起来。

但是又有一个问题——政治问题,地图是给一个国家服务的,那么这地图就要尽可能描述准确这个国家的地形地貌,尽量减小误差,至于别国就无所谓。

所以,就可以人为的把地球的质心“移走”,将局部的表面“贴到”该国的国土,使之高程误差尽量减小到最小。

这个时候,就出现了所谓的“参心坐标系”。即椭球中心不在地球质心的坐标系。如下图:

绿色的球就是为了贴合赤道某个地方而产生了平移的参心系(这里只是个例子,而且画的有点夸张)。

我国常用的参心系及对应椭球:

  • 北京54坐标系:克拉索夫斯基椭球体
  • 西安80坐标系:IAG75椭球体

我国常用的地心系及对应椭球:

  • WGS84坐标系:WGS84椭球体(GPS星历的坐标系,全球统一使用,最新版于2002年修正)
  • CGCS2000坐标系:CGCS2000椭球体(事实上,CGCS2000椭球和WGS84椭球极为相似,偏差仅有0.11mm,完全可以兼容使用)

为什么CGCS2000和WGS84要略微有些偏差?这是因为WGS84系是GPS的坐标系,而我国北斗定位则是需要自己的坐标系,就搞了一波CGCS2000。

这几个坐标系的介绍放在下一节,而这些椭球体的转换将在第三部分介绍(主要就是数学中,空间直角坐标系旋转的问题)。

1.3 我国常见GCS

借助以下4个常见坐标系及椭球体,就可以推及到世界各地不同的GCS及椭球体,完成数据的转化问题。

1.3.1 北京54坐标系(参心)

新中国成立以后,我国大地测量进入了全面发展时期,再全国范围内开展了正规的,全面的大地测量和测图工作,迫切需要建立一个参心大地坐标系。由于当时的“一边倒”政治趋向,故我国采用了前苏联的克拉索夫斯基椭球参数,并与前苏联1942年坐标系进行联测,通过计算建立了我国大地坐标系,定名为1954年北京坐标系。因此,1954年北京坐标系可以认为是前苏联1942年坐标系的延伸。它的原点不在北京而是在前苏联的普尔科沃。

  • 椭球体:Krasovsky椭球
  • 极半径b=6 356 863.0187730473 m
  • 赤道半径a=6 378 245m
  • 扁率=1/298.3
  • 高程系:56黄海系

1.3.2 西安80坐标系(参心)

改革开放啦,国家商量要搞一个更符合国用的坐标系——西安80坐标系,该坐标系的大地原点设在我国中部的陕西省泾阳县永乐镇,位于西安市西北方向约60公里。

  • 椭球体:IAG椭球(全名是啥还得去翻翻课本。。。)
  • 极半径b=6 356 755m
  • 赤道半径a=6 378 140m
  • 扁率=1/298.25722101
  • 高程系:85黄海系

1.3.3 WGS84坐标系(地心)

全称World Geodetic System - 1984,是为了解决GPS定位而产生的全球统一的一个坐标系。

  • 椭球体:WGS84椭球
  • 极半径b=6 356 752.314 245 179 5m
  • 赤道半径a=6 378 137 m
  • 扁率=1/298.257223563
  • 高程系:?根据国家需求定?

1.3.4 CGCS2000坐标系(地心)

2000国家大地坐标系是全球地心坐标系在我国的具体体现,其全称为China Geodetic Coordinate System 2000,其原点为包括海洋和大气的整个地球的质量中心。

  • 椭球体:CGCS2000坐标系
  • 极半径b=6 356 752.314 140 355 8m
  • 赤道半径a=6 378 137m
  • 扁率=1/298.257222101
  • 高程系:85黄海系

【注】CGCS2000的定义与WGS84实质一样。采用的参考椭球非常接近。扁率差异引起椭球面上的纬度和高度变化最大达0.1mm。当前测量精度范围内,可以忽略这点差异。可以说两者相容至cm级水平

最后一张表总结一下:

有趣的是,在ArcGIS的GCS文件夹下,找到了一个“新北京54坐标系”,这是为了使54和80之间方便转化而产生的一个过渡坐标系。

2. 平面坐标与PCS

说完了以经纬度为计量单位的GCS,那么我再来说说以平面(空间)直角坐标系为度量衡的投影坐标系(PCS,Projection Coordinate System)。

说一个具体的问题以解释为什么要用PCS。

如何用经纬度表达一块地的面积?

这没办法吧?经纬度本身不带单位,度分秒仅仅是一个进制。

而且同样是1度经度,在不同的纬度时代表的弧段长是不一样的。

这就给一些地理问题带来了困惑:如何建立一个新的坐标系使得地图分析、空间分析得以定量计算?

PCS——投影坐标系就诞生了。

我要着重介绍一下我国的6种常用投影方式:

  • 高斯克吕格(Gauss Kruger)投影=横轴墨卡托(Transverse Mercator)投影
  • 墨卡托(Mercator)投影
  • 通用横轴墨卡托(UTM)投影
  • Lambert投影
  • Albers投影
  • Web Mercator(网络墨卡托)投影

很多课本、博客都写的很详细了,我想从3D的图形来描述一下他们是怎么个投影的。

2.1 从投影说起

如上图。光线打到物体上,使得物体产生的阴影形状,就叫它的投影。这个不难理解。

这里我想问一个问题:既然投影物体,是不变的,那么我把投影的平面改为曲面呢?

这就产生了不同的投影,比如投射到一个圆锥面上,一个圆柱面上,一个平面上...等等。

不同的投影方式有不同的用途,也有了不同的投影名称。

但是,PCS是基于存在的GCS的,这个直接规定。没有GCS,就无从谈PCS,PCS是GCS上的地物投射到具体投影面的一种结果。

即:

PCS=GCS+投影方式

2.2 我国常见投影

2.2.1 高斯克吕格投影/横轴墨卡托投影

英文名Gauss Kruger。在一些奇奇怪怪的原因中,又名横轴墨卡托投影,英文名Transverse Mercator

它的投影面是椭圆柱面,假设椭圆柱躺着,和地轴垂直,而且投影面与之相切,就是横轴墨卡托了。

中央那条黑线就是投影中心线,与椭圆柱面相切。这条线逢360°的因数就可以取,一般多用3度带、6度带。

就是说,这个投影椭圆柱面可以继续绕着地轴继续转,图中还有一条经线,两条相差6度。

椭圆柱面旋转6度,继续投影,直到360/6=60个投影带投影完毕。

注意3度带和6度带的起算经线不同,以及Y方向(赤道方向)前需要加投影带号。

高斯克吕格已经广为熟知了,我就不作具体介绍,大家可以找比我解释的更好的,我只是摆个图希望大家看的更仔细。

这个投影的特点是,等角/横/切椭圆柱/投影

  • 投影后的地图,角度不变,面积有变。离中央经线越远的地区,面积变化越大。此投影合适用于导航。
  • 投影椭圆柱面是横着的;
  • 投影椭圆柱面与椭球体相切。

适用比例尺:1:2.5万~1:100万等使用6度分带法;1:5000~1:10000使用3度分带法。

【注】在ArcGIS中,不同的GCS的PCS是不同的,以CGCS2000、西安80和北京54为例:

CGCS2000_3_Degree_GK_CM_111E:CGCS2000的GCS下,使用高斯克吕格3度分带法,以中央经线为东经111度的投影带的投影坐标系

CGCS2000_3_Degree_GK_Zone_30:CGCS2000的GCS下,使用高斯克吕格3度分带法,第30个投影带的投影坐标系

Beijing_1954_3_Degree_GK_CM_111E:北京54的GCS下,使用高斯克吕格3度分带法,以中央经线为东经111度的投影带的投影坐标系

Beijing_1954_3_Degree_GK_Zone_35:北京54的GCS下,使用高斯克吕格3度分带法,第35个投影带的投影坐标系

Xian_1980_3_Degree_GK_CM_111E:西安80的GCS下,使用高斯克吕格3度分带法,以中央经线为东经111度的投影带的投影坐标系

Xian_1980_3_Degree_GK_Zone_34:西安80的GCS下,使用高斯克吕格3度分带法,第34个投影带的投影坐标系

不难发现,都是以GCS起头的命名法。

2.2.2 墨卡托投影

英文名Mercator投影。

数学上,投影面是一个椭圆柱面,并且与地轴(地球自转轴)方向一致,故名:“正轴等角切/割圆柱投影”。

既可以切圆柱,也可以割圆柱。

其实就是高斯克吕格的圆柱面竖起来

2.2.3 通用横轴墨卡托投影(UTM投影)

英文全称Universal Transverse Mercator。是一种“横轴等角割圆柱投影

和高斯克吕格类似,高斯克吕格的投影面是与椭球面相切的,这货与椭球面相割。

实质上

  • UTM投影=0.9996高斯投影

其余性质都和高斯克吕格投影一样。

割于纬度80°S和84°N。中央经线投影后,是原长度的0.9996倍。

不过,起算投影带是180°经线,174°W则是第二个投影带的起算经线。

由于有以上优点,UTM投影被许多国家和地区采用,作为大地测量和地形测量的投影基础。

【注】UTM投影是我国各种遥感影像的常用投影。

【注2】UTM投影在ArcGIS中的定义

例如:

WGS_1984_UTM_Zone_50N,就代表WGS1984的GCS下,进行UTM投影,投影带是50N.

WGS_1984_Complex_UTM_Zone_25N,就代表WGS1984的GCS下,进行3度分带UTM投影,投影带是25N.

2.2.4 Lambert投影

中文名兰伯特投影、兰博特投影。

我国地形图常用投影,比如1:400万基础数据:

(GCS是北京54)可以看到授权是自定义,说明这个投影是自定义的,没有被官方收录。等到第三部分再说怎么自定义投影。

我国的基本比例尺地形图(1:5千,1:1万,1:2.5万,1:5万,1:10万,1:25万,1:50万,1:100万)中,1:100万地形图、大部分省区图以及大多数这一比例尺的地图多采用Lambert投影。

兰伯特投影是一种“等角圆锥投影”。

ArcGIS中的投影系一般带有Lambert_Conformal_Conic等字样,国际上用此投影编制1∶100万地形图和航空图。

它就像是一个漏斗罩在乒乓球上:

更标准的画法,见下图,有切和割两种。

它没有角度变形。

这个漏斗的倾斜程度,就有三种:正轴、横轴、斜轴。就是圆锥的方向和地轴的方向的问题。

2.2.5 Albers投影

中文名阿伯斯投影。又名“正轴等积割圆锥投影”,常用于我国各省市的投影。

和上一个兰伯特图形类似,就是一个圆锥与椭球面切割,进行等积投影。

给了官方WKID:102025.

与Lambert投影的区别大概就在一个等角,一个等积投影了。

2.2.6 Web墨卡托(WebMercator投影)

这是一个由Google提出的、为了自家GoogleMap而专门定义的一种投影,是墨卡托投影的一种变种。

主要是将地球椭球面当作正球面来投影,这就会导致一定的误差。

直接看看ArcGIS中的定义:

给了WKID:3857,名字是WGS_1984_Web_Mercator_Auxiliary_Sphere,意思就是在WGS84的GCS下进行web墨卡托投影。

现在,经常被百度地图等网络地图采用,估计是Web程序员想省事吧。

3. GCS与PCS的转换问题(ArcGIS实现)

3.1 GCS转GCS

这就是属于空间解析几何里的空间直角坐标系的移动、转换问题,还有个更高级的说法——仿射变换。

我们知道,空间直角坐标系发生旋转移动缩放,在线性代数里再常见不过了。在摄影测量学中,旋转矩阵就是连接像空间辅助坐标系与像空间坐标系的转换参数(好像不是这俩坐标系,忘了)

欲将一个空间直角坐标系仿射到另一个坐标系的转换,需要进行平移、旋转、缩放三步,可以无序进行。

而平移、旋转又有三个方向上的量,即平移向量=(dx,dy,dz)和旋转角度(A,B,C),加上缩放比例s,完成一个不同的坐标系转换,就需要7参数

我们知道,地心坐标系是唯一的,即原点唯一,就说明平移向量是0向量,如果缩放比例是1,那么旋转角度(A,B,C)就是唯一的仿射参数,即3参数

上图左图为坐标系平移,右图为坐标系旋转。缩放可以在任意阶段进行。

——————以上为理论预备——————

说了这么多理论,如何进行GCS转换呢?假设一个数据源已经有了GCS,我们需要做的操作只有一个:

打开如下工具:数据管理工具/投影和变换/投影,设置界面如下(以WGS1984转西安80为例):

别选错了,这里输入输出都是GCS。然后出现以下警告:

这就告诉你,需要参数转换。在这里,WGS84转西安80,是属于7参数转换(地心转参心),但是缺少7参数,就需要自己去测绘局买或者自己粗略算。

那么如何定义一个地理坐标变换呢?

使用投影的旁边的工具:

即可。见下图:

使用Position_Vector方法(即7参数法)即可输入7参数。我就不输入了,各位有数据的可以继续做。

关于3参数和7参数,在ArcGIS帮助文档里都写有的,目录如下:

3.2 GCS进行投影

这个就更简单了。

随便挑个GCS,喜欢什么用什么,如西安80投影到UTM投影,都可以的。

仍然是上节提及的“投影工具”:

这样就可以了,这里是以WGS84的GCS投影到UTM的第50分度带上。

如果是进行栅格数据的投影,就用“栅格”文件夹下的“投影栅格”工具。

如果所需投影系没有自己需要的GCS,就自定义一个:

这个窗口在Catalog浮动窗或者Catalog软件里打开某个数据的属性,找到XY坐标系的选项卡,就可以新建。

【注】如果在数据的属性页的XY坐标系选项卡,或者图层数据框的XY坐标系选项卡中修改GCS,这仅仅是改个名,坐标值还是原来的坐标系上的,这代表老坐标值并没有转换到新坐标系上。形象的说,就是换汤不换药,这是不对的。我这里说的用投影的方法,才是真正的坐标仿射变换到新的坐标系,使之更改数值,形成在新的坐标系下的新坐标值。

3.3 PCS转PCS(重投影)

最常见的就是下载了谷歌影像图,是Web墨卡托的投影,但是实际又需要高斯投影,那么基于WGS84这个GCS,就可以进行重投影。

在这里,我就以UTM投影转Web墨卡托投影为例:

这次是用“栅格”文件夹下的“投影栅格”工具:

一般选好红框的三个参数即可。

如果仍然提示需要地理坐标变换的警告,说明不是一个GCS的数据,需要3参数或者7参数转换。

栅格数据类似,使用“投影工具”。

工具定位。

3.4 定义投影

这不是定义一个投影坐标系,而是给有坐标值的矢量或者栅格数据添加一个投影坐标系而已。

使用“定义投影”工具即可,既可以定义GCS,也可以定义PCS(这软件的中文翻译有点毛病)。

3.5 地理配准与空间校正

这个就不多说了,地理配准就是使屏幕坐标系的扫描地图仿射、二次三次变换到真正投影坐标系的过程,自动加上目标数据的PCS。有了PCS后就会自动加上GCS。

地理配准主要是针对栅格数据。

空间校正则是针对矢量数据进行仿射、二次、三次变换。

3.6 可能出现的错误

3.6.1 显示几十万位数字的“经纬度”

如上图。

这是有了PCS后,在Catalog的数据属性页的XY坐标系选项卡里,选中GCS,然后应用的结果。

原本是方里网的数字,变成了GCS才有的度分秒。

解决方法:Catalog属性页将GCS改回原来的PCS即可。

3.6.2 显示三位数、两位数的“米”

这个暂时没找到案例,曾经见过。

3.6.3 显示一个几乎是0,一个又很大很大位数的数字

如上图。

这个属于数据本身有GCS,但是在Catalog的XY坐标系选项卡里给它添加PCS然后应用后,可能会出现的错误。

解决方法:在Catalog属性页的XY坐标系选项卡里,选中原来的GCS然后应用即可。

如果数据本身没有PCS,应该做的是投影操作。

3.6.4 大范围的数据给了小范围的投影

例如,整个中国地图理应跨越好几个投影带,却给了某一个投影带的投影坐标系,这就会出现负值。如下图,红框箭头是鼠标的位置。

这个按理说应该用兰伯特投影,但是却给了一个UTM第49区的PCS,所以在中央经线靠左很多的位置会出现负值。

解决方法:这个直接做重投影即可。

 

以上四种错误比较常见,但是手头没有案例,以后遇到再发上来吧。

总结一下:

4. 火星坐标

火星坐标这个东西很常见,出现在互联网地图上。例如百度、腾讯、谷歌等地图。

出于保密等政治因素,地图的GCS坐标值,会被一种特殊的数学函数加密一次,会偏离真实坐标数百米的距离,但是反馈到用户端的却是正确的位置信息(也就是说你拿到GCS坐标也没用,拿GPS到实地跑跟拿着地图定位,可能会偏出几十米甚至一百米的距离)。

火星坐标系原名国测局坐标系(GCJ-02),有篇文章比我写的透彻多了,甚至给出了还原代码,我放到参考资料了,有兴趣的可以看看。

 

参考资料

你所不知的有趣投影方法

ArcGIS将WGS84坐标转北京54或西安80

百度坐标(BD09)、国测局坐标(火星坐标,GCJ02)、和WGS84坐标系互转

&ArcGIS 帮助文档

posted @ 2017-08-30 02:34  岭南灯火  阅读(55356)  评论(33编辑  收藏  举报