地理空间数据的核心就是拥有地理位置信息,而地理位置信息最基本的一个载体就是建立在某一特定坐标系统下的空间坐标。在GIS软件中,空间坐标依据坐标系的不同分为地理坐标和投影坐标。地理坐标是将地球比作一个类椭球体,描述一个点在球面上的位置。但是在地图制图过程中,往往需要在一个平面(无论是纸质地图还是电子地图)上展示地物,这时需要解决地球球面与地图平面之间的矛盾,因此需要对地球进行投影,经过投影后的坐标称为投影坐标,因此投影坐标是建立在地理坐标之上的。
1. 地球椭球体与地理坐标
为了构建地理坐标,首先要理解以下三个基本概念。
(1)地球自然表面。地球自然表面是指由地球上的陆地和海洋表面所构成的自然表面,其形状高低不平,难以构建坐标系。
(2)大地水准面。在假定海水静止的状态下,由海水水面向陆地进行延伸所构成的连续水平面,称为大地水准面。但是,地球质量分布的不均匀性使大地水准面存在轻微的高低不平现象,因此其仍然不适合构建坐标系。
(3)地球椭球面。根据实际需要,构建一个完美的椭球体拟合大地水准面所构成的球面,称为地球椭球面。
地球椭球面是理想的、完美的,地理坐标系也是建立在地球椭球面之上的。但是,地球椭球面不可能完美地拟合大地水准面,因此各个国家或地区建立了能够基本符合自己国家或地区的地球椭球面,或者根据精度需要及特定应用场景构建了不同的地球椭球面。根据构建的地球椭球面的参数不同,地理坐标系也层出不穷。我国的地理坐标系经历了从北京1954坐标系(BJZ54)到西安1980坐标系(XI'AN-80),再到2000国家大地坐标系(CGCS2000)的发展过程。
地理坐标系可以通过WKT语句描述。例如,利用WKT语句描述WGS 1984坐标系:
“GEOGCS”声明一个地理坐标系,第1项为地理坐标系的名称;第2项“DATUM”描述一个坐标基准,并通过“SPHEROID”描述其名称及地球椭球体的基本参数;第3项“PRIMEM”标识一个起始经线,此处以格林尼治0°线为起始经线;第4项“UNIT”定义了地理坐标系的单位;第5项“AUTHORITY”标识地理坐标系编码,这个地理坐标系可以用“EPSG:4326”编码表示。
2. 地图投影与投影坐标
为了解决地球椭球面和地图平面之间的矛盾,需要将地球椭球面进行投影,经过投影以后的坐标系称为投影坐标系。投影后的平面坐标系一定会出现变形,我们只能在等距、等积和等角之间进行取舍。因此,在不同应用场景下,大量的投影坐标系应运而生。
WKT语句也可以用于描述投影坐标系,例如,利用WKT语句描述WGS 1984 Web墨卡托投影坐标系(辅助球面):
“PROJCS”声明一个投影坐标系,第1项为投影坐标系的名称;第2项“GEOGCS”指明该投影坐标系使用的地理坐标系;第3项“PROJECTION”描述投影方法;第4~8项“PARAMETER”声明各种投影参数;第9项“UNIT”定义坐标单位;第10项“AUTHORITY”标识投影坐标系编码,这个投影坐标系可以用“EPSG:3857”表示。
【小提示】 WGS 1984 Web墨卡托投影坐标系(EPSG:3785)与WGS 1984 Web墨卡托投影坐标系(辅助球面)(EPSG:3857)非常类似,后者采用辅助球面将椭球体近似转换为一个正球体后进行投影变换。前者覆盖整个地球范围,后者只覆盖北纬85°与南纬85°之间的范围。
另外,Web墨卡托投影并不是真正的墨卡托投影,而是伪墨卡托投影(Popular Visualization Pseudo Mercator)。
3. EPSG
由于地理坐标系和投影坐标系众多,如果仅通过参数对这些坐标系进行整理与应用则过于麻烦,因此需要通过标准化组织将这些坐标系归档整理。对于石油的探查和开采来讲,坐标系的不同会显著影响开采精度,因此欧洲石油调查组织(European Petroleum Survey Group,EPSG)整合了绝大多数常用的坐标系,并为每个坐标系设置了一个编码,例如,“EPSG:4326”和“EPSG:3785”分别表示WGS 1984坐标系和WGS 1984 Web墨卡托投影坐标系。常用的坐标系及其EPSG编码如表1-1所示。
表1-1 常用的坐标系及其EPSG编码
各坐标系的EPSG编码可在这一网站查找:http://spatialreference.org/ref/epsg/。
4. GCJ-02坐标系
GCJ-02坐标系是由中国国家测绘局制定的中国范围内民用地图(包括电子地图)加密后的地理坐标系,字母“GCJ”分别为“国家”、“测绘”和“局”的首字母简写。GCJ-02的加密算法是非线性的,很难通过加密坐标反推原始的坐标。目前,谷歌、百度等公司发布的电子地图基本采用GCJ-02坐标系对原始数据进行加密。因此,通过GPS获取的WGS 1984坐标系的位置信息经过投影后无法直接准确地叠加在电子地图上,需要经过GCJ-02加密处理。
地理要素经过抽象才能存储在计算机中,按照存储数据结构的不同,其可分为矢量数据和栅格数据两类(见图1-2)。
图1-2 矢量数据结构和栅格数据结构
· 矢量数据利用记录地理要素坐标(也可包括拓扑关系)的方式进行存储,并将属性信息放置在单独的属性表中,属性表中的每条记录都与地理要素一一对应,具有定位明显、属性隐含的特点。
· 栅格数据将一个空间划分成一系列规律分布的格网,并用每个格网的数值表达属性信息,这些格网构成的数据阵列通过参数设置的方式将其放置或拟合在特定的坐标系中,具有属性明显、定位隐含的特点。
相对来说,矢量数据的空间位置精度比较高,并且输出的地图一般比较容易且细腻,但是由于矢量数据的录入与管理比较复杂,有时需要维护其拓扑信息,因此往往给数据管理和空间分析带来一些麻烦。栅格数据一般在空间分析操作上较简单,并且不需要维护其拓扑关系,但是常常占据较大的存储空间,而且空间位置精度较低。
网格数据(Mesh Data)是指在二维空间或三维空间中,通过顶点(Vertices)、边(Edges)与面(Faces)等方式记录多维数据,是矢量数据和栅格数据的补充形式,多用于存储气候气象数据、水文数据、洋流数据等。
1. 矢量数据
矢量数据通常用点、线、多边形等基本几何对象描述地理要素。几何对象可以通过WKT和WKB两种方式进行描述。
1)WKT
WKT采用文本的方式存储几何对象,主要包括类型声明部分和数据部分。其中,类型声明部分采用英文文本进行声明,数据部分放置在类型声明部分的后面,并用小括号括起来。对于具有多个节点的几何对象(如线、多边形等),同一个节点内的数据用空格隔开,不同节点之间的数据用逗号隔开,例如,一个坐标为(20,15)的点对象的声明方式如下:
“Point”表示该几何对象的类型为点对象,小括号中的第一个数字表示 X 坐标,第二个数字表示 Y 坐标。
包含Z值和M值的点需要在“Point”之后通过字母“Z”与“M”进行声明,并将Z值和M值加入小括号中,例如:
“Point Z”中的“10”表示Z值,“Point M”中的“8”表示M值,“Point ZM”中的“10”和“8”分别表示Z值和M值。
【小提示】 任何一个几何对象的任何一个节点都可以包括Z值和M值,Z值通常表示高程信息,将几何对象推向三维;M值通常表示其他属性信息,如温度、湿度等,将几何对象推向第四个维度。
利用WKT表达一个具有四个节点的线,其节点坐标分别为(10,10),(20,20),(30,30)和(40,40):
利用WKT表达一个具有四个节点的多边形,其节点坐标分别为(10,10),(10,20),(20,20)和(20,15):
注意,最后一个节点和第一个节点相同,以闭合图形。
在上述WKT代码中,“Polygon”后面由两个小括号组成,第一层小括号标识数据为“Polygon”数据;第二层小括号标识一个环(Ring),包括内环(Inner Ring)和外环(Outer Ring)两类。一个多边形至少包括一个外环。利用内环可以将一个多边形“挖出”中空的部分;利用多个外环可以将一个多边形分为多个部分,表达“飞地”的地理要素概念。
例如,当需要一个多边形表达两个岛屿组成一个地理要素单位,并且其中一个岛屿具有中空的“水体”时,可以利用两个外环和一个内环组成一个多边形,如图1-3所示。
图1-3 由两个外环和一个内环组成的多边形
另外,利用几何对象集合(GeometryCollection)可以将多个不同类型的几何对象组合为一个几何对象。例如,将一个坐标为(10,10)的点和由两个节点组成的线组合为一个几何对象,WKT代码如下:
2)WKB
WKB采用二进制方式存储几何对象信息,包括字节序、几何类型和坐标三部分。· 字节序部分:可以为0(大端字节序,Big-Indian)或1(小端字节序,Little-Indian),占一个字节。
· 几何类型部分:用4个字节的编码表示一个几何类型,几何类型和编码的对应关系如表1-2所示。
· 坐标部分:一个坐标( X 坐标、 Y 坐标、Z值、M值)用一个双浮点类型(Double)空间存储,占8个字节。
表1-2 WKB几何类型和编码的对应关系
例如,采用WKB方式表示一个坐标为(20,15)的点的十六进制:
这个点对象共占据21个字节,第1个字节表示大端字节序;第2~5个字节“00 00 00 01”表示几何类型为点;第6~13个字节“40 34 00 00 00 00 00 00”表示 X 坐标20,第14~21个字节“40 2E 00 00 00 00 00 00”表示 Y 坐标15。
可见,采用WKT方式描述几何图形更直观,但是由于WKT字符串需要被解释器解释后才能被计算机处理,所以读取速度相对较慢;采用WKB方式描述几何图形更适合计算机读取,但是失去了直观性。
2. 栅格数据
栅格数据采用某种数据类型的数值阵列存储数据,阵列中的每个数值称为一个像元(Pixel)。由于数据阵列本身不存在空间信息,因此需要元数据进行界定。栅格数据的元数据包括空间坐标系、数据类型等。
在ENVI DAT栅格数据格式中,后缀名为“hdr”的头文件存储了完整的元数据信息。此处以某ENVI DAT格式数据的“hdr”文件为例,介绍常用的元数据信息:
1)数据类型与行列数波段数
数据类型是指一个栅格像元中数值的数据类型(Data Type)。GDAL规定了12种栅格数据类型,可以覆盖绝大多数的栅格数据(见表1-3)。
表1-3 常用的栅格数据类型及其GDAL定义
栅格数据的行列数是指栅格像元阵列的行列数。注意,栅格数据行列号的定义是从左上角开始的(栅格空间),因此第一行位于栅格像元阵列的最上方,第一列位于栅格像元阵列的最左方,这与坐标系的方向定义(坐标空间)不同(见图1-4)。
图1-4 栅格空间行列号的定义与坐标空间的方向定义
2)坐标参考系
坐标参考系(Coordinate Reference System,CRS)界定了栅格数据所处的投影坐标系或地理坐标系。
3)参考坐标与像元大小
GeoTiff、PNG等格式的栅格数据常常存在一个单独的世界文件(World File),它通过六个参数界定栅格影像的行列号与实际地理坐标的仿射关系,通常称为栅格数据的“六参数”。
世界文件的名称通常与栅格数据文件的名称相同,但其后缀名是在栅格数据文件的后缀名后加上字母“w”,如“jpgw”、“tifw”、“tiffw”和“pngw”等。
但是,有时为了保证文件扩展名为三个字符,也可以把世界文件扩展名中间的字母去掉,如“jgw”、“tfw”和“pgw”等。
另外,世界文件的扩展名也可以设置为“wld”。
【小提示】 世界文件并不是必需的,因为这些参数还可以保存在栅格数据的元数据文件(如ENVI DAT格式数据的“hdr”文件)中,也可以保存在数据本身(如GeoTiff文件的内部标签)中。
世界文件一共有六行数值,每行数值代表一个参数,各定义如下:
· 第一行(A):栅格数据中一个像元的宽在 x 轴方向的大小,使用地图坐标单位。
· 第二行(D):栅格数据中一个像元的宽在 y 轴方向的大小,使用地图坐标单位。
· 第三行(B):栅格数据中一个像元的高在 x 轴方向的大小,使用地图坐标单位。
· 第四行(E):栅格数据中一个像元的高在 y 轴方向的大小,使用地图坐标单位。由于栅格空间和坐标空间在纵向方向的相反关系,该参数通常为负值。
· 第五行(C):栅格数据左上角像元中心点的 x 坐标,使用地图坐标单位。
· 第六行(F):栅格数据左上角像元中心点的 y 坐标,使用地图坐标单位。
地图坐标与行列号的仿射关系为:
x = A ×col+ B ×row+ C
y = D ×col+ E ×row+ F
其中, x 和 y 分别表示 X 坐标和 Y 坐标;col和row分别表示列号和行号。
4)坐标范围
栅格数据的坐标范围(Extent)通常由 X 坐标最小值(xmin)、 X 坐标最大值(xmax)、 Y 坐标最小值(ymin)、 Y 坐标最大值(ymax)组成,它们通常被称为栅格数据的“四至范围”。
5)栅格数据的存储格式
栅格数据的存储格式包括波段顺序格式、波段按行交叉、波段按像元交叉等类型。
(1)波段顺序格式(Band Sequential Format,BSQ)采用按波段存储的方式,先存储一个波段中的数据,然后存储另一个波段中的数据,直到全部数据存储完毕。
(2)波段按行交叉(Band Interleaved by Line Format,BIL)采用按行存储的方式,先存储第一行的第一波段数据,然后存储第一行的第二波段数据,直到第一行的数据全部存储完毕;然后按照上述方式存储其他行数据,直到全部数据存储完毕。
(3)波段按像元交叉(Band Interleaved by Pixel Format,BIP)采用按像元存储的方式,先存储第一个像元全部波段的数据,然后存储其他像元全部波段的数据,直到全部像元存储完毕。
上述三种存储格式的读取性能存在差异,BSQ适合单个波段中部分区域数据的存储,BIP适合对某个像元的波谱信息进行存取,BIL则是BSQ和BIP的折中方式。当栅格数据只包括一个波段时,BSQ、BIL和BIP的存储没有区别。
3. 网格数据
作为矢量数据和栅格数据的补充,网格数据兼顾了两者的优势:既可以像矢量数据一样在顶点上存储数据(Defined on Vertices),也可以像栅格数据一样在面上存储数据(Defined on Face)。网格数据既可以像矢量数据一样以不规则的多边形描述数据的定位信息,也可以像栅格数据一样以规则的网格描述数据的定位信息,如图1-5所示。
图1-5 网格数据的网格类型
网格数据存储大量信息的同时占据较少的磁盘空间,例如,存储高程信息的TIN数据通常比栅格存储方式的数据要小很多,并保持较高的精度。网格数据具有特定的应用方向,主要包括:
(1)存储具有空间自相关(满足地理学第一定律)的数据。
(2)存储长时序(包含多个时刻)的数据。
(3)存储具有矢量(方向)的数据,如风向、流向等。
因此,网格数据多应用在气象、水文等领域。
在QGIS中,网格数据的读取依赖网格数据抽象库(Mesh Data Abstraction Library,MDAL)。与GDAL类似,MDAL采用C++语言编写,并通过MIT协议发布。MDAL起源于QGIS 2中的Crayfish插件,目前主要应用在QGIS软件中,并为MDAL提供大量反馈信息。自QGIS 3.2以来就逐渐集成了MDAL,并在QGIS 3.4逐渐稳定成熟。MDAL支持2DM、NetCDF、GRIB等多种网格数据格式,其支持的全部数据类型可以在https://www.mdal.xyz/网站中查询。
【小提示】 在QGIS 3.10中,MDAL仅支持1D和2D网格数据;在QGIS 3.12中,MDAL开始支持3D网格数据。
地理空间数据的尺度可以分为空间尺度和时间尺度,空间尺度是指数据在空间上的精细程度,时间尺度是指数据在时间上的精细程度。在地图制图和空间分析中,应当注重数据的尺度选择,本节介绍地理空间数据的尺度问题及相关的重要概念。
1. 尺度效应与尺度选择
尺度效应是指空间数据分析和展示结果会因数据在空间或时间上的最小信息单元水平的变化而存在差异(见图1-6)。在实际的数据分析中,我们需要根据分析或研究的目的,选择适合这项分析或研究尺度的数据,即尺度选择。数据的尺度不宜过高,例如,当我们使用遥感数据研究城市道路时,需要选择高分辨率遥感数据,如果选择MODIS这类资源遥感卫星数据,那么道路无法有效地被提取,其特征将被淹没在混合像元(Mixed Pixel)中。数据的尺度也不宜过低,例如,当我们进行中国的森林变化这样大尺度的研究时,就不宜选择高分辨率遥感数据,否则会带来大量低效的时间成本和计算成本。无论是矢量数据还是栅格数据,在数据的采集、编辑和处理的过程中,一定要时刻注意尺度效应,这样才能事半功倍。
图1-6 不同尺度的栅格数据
在栅格数据中,空间分辨率不仅意味着精度不同,同种数据在不同尺度反映的地理问题往往也是不同的。
2. 栅格数据金字塔
栅格数据金字塔是在不同尺度创建数据副本,以快速地在某个尺度展现数据,即利用冗余存储的方式换取数据的响应速度。数据金字塔的底部是栅格原始数据,从下至上是分辨率逐渐降低的数据快照。
3. 制图综合
制图综合(Cartographic Generalization)是地图学中的重要概念,是指在从大比例尺地图转换为小比例尺地图的过程中,根据制图的实际需求与目的,将地理要素加以概括和抽象,留下最关键的信息,从而突出制图者希望传达的主要内容。
希望读者在学习中特别关注数据的尺度效应问题,这样才能制作出精美的地图,快速得到想要的数据结果。