陈 静,陈帆飞,周 强
(武汉大学测绘遥感信息工程国家重点实验室,湖北武汉430079)
全球地形数据多尺度晕渲及服务研究
陈 静,陈帆飞,周 强
(武汉大学测绘遥感信息工程国家重点实验室,湖北武汉430079)
采用地形分块晕渲的策略,根据不同地形尺度,采用不同的晕渲方法对全球地形数据进行多尺度晕渲。在此基础上,以WCS方式提供全球地形晕渲数据服务,最后通过试验对上述方法的有效性进行验证。
全球地形数据;分块晕渲;WCS
目前,地形数据的自动晕渲还不能完全满足数字地图制图的要求,主要表现在自动晕渲还不能完全表达所有的地貌形态,不能进行局部光照调整,不能快速处理海量数据等。对此,国内外相关专家做了大量的研究探索工作[1-3]。对于数字晕渲法的已有研究主要集中在两个方面:①对于地形表面的模拟以及法向量的计算;② 对于光照的模拟与计算。
已有研究中只涉及局部地区地形地貌的晕渲,且多用于局部直角坐标系下的地形晕渲,尚没有涉及全球大范围地形晕渲的研究。
由于全球地形数据具有海量的特征,进行全球范围的地形晕渲计算时,若一次性调入内存进行计算则要浪费大量的系统资源,甚至有可能造成系统崩溃。同时,由于全球地形数据具有多尺度特征,不同的尺度条件下,地形地貌的特征不同,而晕渲图又是一种直观的对地形数据进行二维可视化的方法。因此,对多尺度全球地形数据,若采用同一种晕渲方法,则无法突出不同尺度下晕渲表达的重点。最后,生成的全球多尺度晕渲图的海量数据也需要进行高效地组织与管理,并将其封装成网络覆盖服务(Web coverage service,WCS),实现地形晕渲的信息共享。
本文在设计全球地形数据逻辑统一索引的基础上,重点探讨了基于分块策略的多尺度地形晕渲方法,并将晕渲结果封装成WCS,提供地形晕渲服务。
由于全球地形数据具有海量、多尺度、多时相等特征,因此无法对如此巨大的地形数据进行一次性晕渲,需要考虑按照一定范围对全球地形数据进行分块晕渲,并且在逻辑上组织成全球统一的地形晕渲图。对此,设计了全球地形数据多尺度晕渲的逻辑分块索引方法,其核心是面向全球范围,基于Ottoson提出的用椭球四叉树对全球地理数据进行索引的方法[4],构建全球等经纬度的地理坐标系,分别对东西半球建立两个层次四叉树结构。四叉树4个结点分别代表该逻辑块的空间范围,全球四叉树逻辑索引一旦建立,全球地形数据多尺度晕渲的逻辑索引也就确定下来,如图1所示。
图1 全球范围逻辑金字塔结构图
全球地形数据多尺度晕渲的逻辑索引实际上是经纬网点阵列,每一个格网与其相邻格网之间的拓扑关系隐含在该阵列的行列号之中。对此,规定以每一个逻辑分块的左下角行号、列号构成空间索引关键字,将行号和列号编码成Morton码。同时考虑到地形数据的时态信息,扩展上述Morton码与时间编码组成一个分块的时空统一编码,用一个64 bit整型表示。
全球地形多尺度晕渲方法主要包括:①基于分块的地形晕渲;②晕渲的色调统一处理;③多尺度晕渲。
1.基于分块的地形晕渲
在全球多尺度地形数据晕渲的逻辑索引上,对逻辑索引中任意尺度下某一个逻辑块,通过Morton码计算出其经纬度范围,并根据此范围和尺度,调度相应尺度和范围的原始地形数据,进行分块的地形晕渲。为了便于晕渲结果的组织与管理,规定这个逻辑块范围也是最终生成晕渲结果图的范围。
具体来说,大范围分块晕渲,首先是根据地形数据最大分辨率的大小,找出与其最接近的全球索引的层号,然后根据式(1),计算出需要晕渲的范围在全球索引层中的行列号范围。最后根据块的行列号范围,逐层逐块进行晕渲。
式中,Lon和Lat分别为经、纬度;TileSize为索引块的大小;abs()表示取绝对值;int()表示取整。
在晕渲过程中,由于每个高程点的坡度坡向或者平均法向量都需要周围8个相邻点的高程进行计算,但是在一个逻辑分块的边缘只能取到小于等于4个的高程点数。因此两个相邻逻辑块的边缘计算出来的值很难完全相同,由此会导致块与块之间晕渲像素的不连续。于是,在分块晕渲中,需要取稍大于逻辑晕渲块地理范围的高程数据。即如果晕渲块实际地理范围的起始点为(x0,y0),采样点的个数为m·n(m为采样点的列数,n为采样点的行数)。那么在取数据时,应在原数据的基础上扩大地理范围,起始点变成(x0',y0'),采样点个数为(m+2)·(n+2),即在原采样点阵的4个方向上各增加一行和一列,如图2所示。斜线部分表示的是逻辑晕渲块范围,也是最终生成的晕渲图范围,外面扩展框则是实际取出的用来晕渲的地形数据范围。实际进行晕渲计算时,只计算中间的m·n个点(即图2中斜线范围)的晕渲值。
2.分块晕渲的色调统一处理
由于单个地形数据块的晕渲只在一个数据块内进行,所以不需要考虑色调和光照不一致的问题。而对于同一尺度下大规模分块晕渲的地形来说,必须从整体考虑色彩设定与配置。
若仅参照每个晕渲块内的最大最小高程值对,由于各个块之间的最大最小高程不尽相同,那么晕渲后,块与块之间必然会出现色调不一致的情况。因此,为了避免这种情况必须将绝对高程值进行分级赋值处理。
图2 基于分块的晕渲方法
可以按地形数据的绝对高程将地貌形态划分为几级,每一级根据高程继续细分若干子层,根据经验值给地形数据赋色彩值,每一级内部子层的色彩变化趋势不变。结合该种设色方式,若高程值H在高程分级值Hm和Hm+1之间,则具体的晕渲色彩中RGB 3个颜色分量时运用式(2)进行计算
式中,Rm+1、Rm、Gm+1、Gm、Bm+1、Bm分别是第m+1级和m级高程对应的RGB值;k是基本光照模型[5]中关于地表曲面微元光线照度的一次函数,k具体取值可根据经验和晕渲效果进行调整。
色层表中每一级高程阈值的取值,可以根据具体的地形晕渲范围中高程的分布进行取值。即在进行具体的晕渲之前,先对整个范围内的地形高程值抽样做一个统计直方图,根据统计结果来确定划分高程分级的级数以及每一级的阈值,如图3所示。对于阈值的确定:①首先找出直方图中所有极小频数的高程值,即直方图中下凹的区域;② 其次从最小的极小频数高程值开始,取该值为第一个阈值;③再统计每两个相邻极小频数高程值之间高程的频数,若该频数大于一定的频数值,则取较高的极小频数高程值作为阈值,继续向下统计,若该频数小于一定的频数值,则舍弃较大的高程值;④然后把每一级按高程的频数划分为若干层;⑤最后参照以往的经典分层设色值,如Nighbert的分层设色,对分好的高程层进行设色。
图3 高程直方图
经过上述处理过程的分块晕渲,不仅能够保证块之间色彩色调的一致性,而且由于经过高程直方图的统计,晕渲的层次感较强。
3.地形数据的多尺度晕渲
一般来说,计算晕渲的方法有两种:① 通过坡度、坡向进行计算;②通过该点周围地面三角形片的法向量进行计算,即平均法向量计算法。其中,坡度、坡向计算方法对地形地貌的晕渲表现比较简洁干净、重点地貌形态突出,而通过法向量计算的方法对地形地貌的表现比较破碎,但其细部特征保留较好。因此对于小比例尺的地图的晕渲应该使用对整体地形地貌概括性强的坡度、坡向的方法进行计算,而对于大比例尺和需要对地形地貌精确表达的地图,则应该使用对细节保留较好的平均法向量方法进行计算。
为了方便确定比例尺的尺度,对于地形晕渲方法的选择,参照下述统计判定规则
式中,p、k、n和m为常数,分别为0.004 29、0.000 18、1.444605和0.26777;Scale表示国家基本地形比例尺的分母值除以10 000,如1∶1 000 000,则Scale取100;Res则是该比例尺对应的通用分辨率,如1∶1 000 000对应1 000 m的网格大小;u是坡向信息失真因子[6]。当Res<200,u=0.1;Res在200~1 000,u=0.2,Res>1 000,u=0.3。
当Y<0.5时,采用法向量进行计算,否则采用坡度、坡向计算方法。
地形晕渲图是地形数据的一种直观、简便的可视化表达,将其封装成WCS有利于实现地形晕渲信息的共享与服务。本文将实现WCS的3个操作:对GetCapabilities、DescribeCoverage和 GetCo-verage进行了封装,并在内部设计了从晕渲结果数据集到WCS的数据流程。同时设计了3层架构的服务体系,包括数据层、应用服务层和应用层,提供全球地形多尺度的晕渲数据服务,如图4所示。
图4 基于WCS的地形晕渲服务实现模式
图4中,数据服务层是多尺度晕渲图结果数据,它在Web服务器接口层或者GIS应用服务器中进行注册,并提供相应的元数据。应用服务层包括WCS服务器、GIS应用服务器、Web服务器等。WCS服务器负责WCS接口的封装,处理响应用户的请求,向数据层请求数据和返回WCS的地形晕渲图;GIS应用服务器、Web服务器等主要提供地形晕渲数据的注册服务、目录服务、客户端的数据请求等一系列服务。应用层支持应用程序和IE浏览器的应用,用户可以进行相应的操作和请求。
基于上述架构,采用全球90 m分辨率STRM地形数据,进行全球地形数据多尺度分块晕渲,试验效果图如图5所示。封装成WCS,利用IE浏览器请求得到的结果如图6所示。
图5 多尺度全球地形晕渲图
在保证晕渲色调的一致性条件下,多尺度分块地形晕渲的方法能突出不同尺度下地貌的表达。同时将晕渲结果封装成WCS提供全球地形多尺度晕渲服务。
图6 地形晕渲WCS结果
[1] 施祖辉.地貌晕渲法[M].北京:测绘出版社,1985.
[2] 陈艳丽,李少梅.基于坡度坡向的地貌晕渲实现研究[J].测绘科学,2008(S1):181-182.
[3] 郭庆胜,王晓延.地貌晕渲中光源使用方法与用色规则的研究[J].武汉大学学报:信息科学版,2004,29(1):20-23;33.
[4] OTTOSON P,HAUSKA H.Ellipsoidal Quadtrees for Indexing of Global Geographical Data[J].Geographical Information Science,2002,16(3):223-226.
[5] 韩明峰,Phong光照模型中单位反射光线矢量的快速计算研究[J],微机发展,1999,9(1):6-7.
[6] 刘春,孙伟伟,吴杭彬.DEM地形复杂因子的确定及与地形描述精度的关系[J].武汉大学学报:信息科学版,2009,34(9):1014-1019.
Study on Method and Service of Multi-scales Relief Shading for Global Terrain Data
CHEN Jing,CHEN Fanfei,ZHOU Qiang
0494-0911(2011)07-0027-03
P283.7
B
2010-08-04
国家863计划资助项目(2009AA12Z229);国家863计划资助项目(2004AA12Z201)
陈 静(1975—),男,江西樟树人,副教授,主要从事基于网络环境的地球空间信息服务研究与应用。