基于图像集的大规模数字化高程数据模型自动生成方法.pdf

本发明公开一种基于图像集的大规模数字化高程数据模型自动生成方法，首先将各个数字位图图像缩放为标准化图像，并对其进行半色调处理；再根据输入参数对平面点集进行简化处理，并为简化点集中的每个点设定一个随机高度值及一个圆形控制区域，并在每个控制区域上定义一个控制函数；在每块标准化图像上构造地形网格，并确定网格点的控制点；用控制函数确定的权值将这些控制点的高度值混合在一起，作为虚拟地形块网格点的高度值；最后通过拼合相邻地形块DEM数据产生大规模DEM数据集。本发明实现了不同分辨率的大规模数字化DEM的自动化合成，其原理简单、易于实现，并提供了多种参数对地形合成过程进行控制，能够保证相邻地形块间的自然顺滑过渡。

1.  基于图像集的大规模数字化高程数据模型自动生成方法，其步骤是：
a、对系列位图图像进行标准化预处理：将所有用于生成数字化地形的位图图像都缩放为长宽一致的标准化图像；
b、对标准化图像进行半色调处理：将具有连续色调的标准化图像转化为图像区域内由黑色点构成的、疏密度随空间位置变化而变化的平面点集；
c、简化平面点集：根据设定的点集简化控制参数，运用点集简化方法对平面点集进行简化处理，将平面点集简化为指定密度的简化点集，并为简化点集中的每个控制点赋予一个随机高度值，使该高度值的取值在输入参数所确定的范围内；
d、分解标准化图像并确定子块间的相邻关系：将相关图像区域拼接在一起，使之在平面上铺嵌出一块连续区域；接着，将每块图像区域一分为四，得到4块大小相等、形状相同的子块(LU、LD、RU、RD)；每个子块的相邻子块共有8个，其中3个位于其自身所属的图像区域内，其它5个位于其所属图像区域的相邻图像区域内；
e、确定控制点的控制域和控制函数：控制点(v_i)的控制域是以控制点(v_i)为圆心、以控制域半径(R_i)为半径的圆形区域，其中控制域半径(R_i)由输入的参数确定；简化点集中各控制点的控制域半径(R_i)的最大值不超过图像区域边长的一半；控制函数是一个定义在控制域上的径向单调减函数，其在圆心处取到最大值，随着与圆心距离的不断增加，其值不断衰减，并在到达控制域边界时衰减为0；
f、构造地形栅格：在每块图像区域内产生正交地形网格，使该网格恰好覆盖整个图像区域；以用户输入的参数确定网格的密度，再根据网格点在图像区域内的位置，确定其坐标；不同图像区域上的网格具有相同的结构和密度。
g、确定地形网格点的高程值：对于任一网格点(p)，首先确定它所属于的图像区域子块(B₀)，从图像区域子块(B₀)的9个相邻子块(B₀，B₁，B₂，…，B₈)中寻找控制域覆盖网格点(p)的所有控制点(v_i(i＝1，2，…，k))，然后根据控制点(v_i)的高程值(h_i)、控制点(v_i)的权(w_i(p))、以及用户输入的全局控制参数(η)计算网格点(p)的高程值
H(p)=η1kΣi=1kwi(p)hi]]>
计算出所有网格点的高程值后，即得到相应图像区域的地形块DEM；
h、建立大规模地形的DEM模型：分别为每个地形块DEM各设定一个偏移坐标，用以确定该地形块DEM在全局地形坐标系中的位置；各偏移坐标分别为相应图像区域最左下网格点的全局坐标；将地形块DEM的偏移坐标与块内各网格点的局部坐标相加，便可得到各网格点的全局坐标，并让相邻地形块DEM的边界网格点重合，便可构成大规模DEM地形模型数据集。

基于图像集的大规模数字化高程数据模型自动生成方法
技术领域
本发明公开一种利用系列半色调图像自动生成大规模数字化虚拟地形高程数据模型DEM(Digital Elevation Model)的方法，内容涉及图像半色调处理与数据加权混合技术，尤其涉及一种数字化虚拟地形的自动制备方法。
背景技术
数字化地形通常以高程数据模型DEM(Digital Elevation Model)作为表示形式，现有技术公开了多种数字化地形的生成方法。
一类通过测量技术直接由现实世界获得相关地域数字化DEM的生成技术。该类技术由地形测量数据，如由遥测卫星采集到的GIS(Geographic InformationSystem)数据等，直接生成地形的数字化DEM。这类方法的不足是，获取符合特定需要的地形模型相对困难，且获取成本较高，难以适应诸如虚拟驾驶训练以及电子游戏等应用的需求。
一种称为Sculpting的数字化地形生成技术。该技术通过手工塑造微型的物理地形模型，再用三维采集设备对其进行数字化处理来产生所需的DEM。使用该技术构造地形DEM的时间周期相对较长。
一类基于分形过程的数字化虚拟地形DEM生成技术(参见：Lewis J P.Generalized Stochastic Subdivision.ACM Transactions on Graphics，1987，6(3)：167-190)。该类技术通过迭代的随机分形过程，产生具有多尺度自相似性的分形地形外观。该类技术的不足是，很难对地形的生成过程进行有效的控制。后来，人机交互手段被引入到此类技术中，使人们能够通过鼠标等交互设备对地形生成过程进行编辑和修改，但对于大规模数字化地形的生成过程来说，人机交互工作过于繁重。
一种应用于三维游戏的大型地形生成方法(参见：王建民等.专利“一种基于PC平台的三维游戏中大型地形生成方法”，申请号：CN200610122654.7)和一种基于视点的地形网格生成方法(参见：华炜等.专利“一种与视点相关的基于四队列优化的地形网格生成方法”，专利号：ZL200310108882.5)。这两种方法都提出了如何由地形文件数据来生成用于显示的大型数字化虚拟地形模型，但二者都需要用已经制备好的地形DEM文件作为方法的输入，其方法的核心在于解决地形绘制过程中的网格动态管理问题，而非解决地形DEM模型的合成与制备问题。
一种基于图像径向权混合的大规模虚拟地形合成技术(庞明勇等.大规模虚拟地形的图像径向权混合生成算法。南理工大学报(自然科学版)，2009，33(4)：438-443)。该技术以一组图像作为输入，在对图像进行滤波去噪等处理的基础上，将图像像素的灰度值作为地形高程值，再通过有限支撑函数的混合操作把不同图像沿边界处融合，产生大规模数字化虚拟地形。该技术的不足之处在于，由于混合操作中对地形子块进行了对称翻转叠加等处理，使相邻地形块在拼接处附近产生了对称性明显的特征，从而在生成的地形中留下人工加工的痕迹。
发明内容
本发明的目的是针对现有技术的不足，提出一种基于图像集的大规模数字化高程数据模型自动生成方法，是一种基于图像半色调处理与数据加权混合技术、利用图像信息自动生成DEM的方法。生成的大规模DEM由一块块地形块拼合而成，各个地形块可以独立存储，相邻地形块之间的过渡自然平顺，且无人工加工痕迹。
为实现上述目的，本发明采用以下技术方案：
1、对系列位图图像进行标准化预处理：对尺寸大小不同的数字位图图像进行缩放处理，将所有用于生成数字化地形的位图图像都缩放为长宽一致的正方形标准图像(称为标准化图像)。标准图像占据的平面区域称为图像区域。
2、对标准化图像进行半色调处理：采用半色调处理技术将标准化图像转化为黑白二值点表示的半色调图像，即将具有连续色调的标准化图像转化为图像区域内由黑色点构成的、疏密度随空间位置变化而变化的平面点集。
3、简化平面点集：根据设定的点集简化控制参数，运用点集简化方法对平面点集进行简化处理，将平面点集简化为指定密度的新点集(称为简化点集)，并为简化点集中的每个点(称为控制点)赋予一个随机高度值，使该高度值的取值在输入参数所确定的范围内。
4、分解标准化图像并确定子块间的相邻关系：将相关图像区域拼接在一起，使之在平面上铺嵌出一块连续区域。接着，将每块图像区域一分为四，得到4块大小相等、形状相同的子块，分别记为LU、LD、RU、RD(参见附图14)。每个子块的相邻子块共有8个，其中3个位于其自身所属的图像区域内，其它5个位于其所属图像区域的相邻图像区域内。如附图15所示，共有A、B、C、D四个相邻的图像区域，其中A  中子块A_RD的8个相邻子块分别为A_LD，A_LU，A_RU，B_LU，B_LD，D_LU，C_RU和C_LU。对于其它子块，类似地确定它们的相邻子块。
5、确定控制点的控制域和控制函数：控制点v_i的控制域是以控制点v_i为圆心、以控制域半径R_i为半径的圆形区域，其中R_i由输入的参数确定。简化点集中各控制点的控制域半径R_i可以相等，也可以不同，但其最大值不能超过图像区域边长的一半。控制函数是一个定义在控制域上的径向单调减函数，其在圆心处取到最大值，随着与圆心距离的不断增加，其值不断衰减，并在到达控制域边界时衰减为0。附图10中，半径为R的圆为控制点v的控制域；附图11示出了一个控制函数的图像。
6、构造地形栅格：在每块图像区域内产生正交地形网格，使该网格恰好覆盖整个图像区域；以用户输入的参数确定网格的密度(即网格中小正边形的数目)。再根据网格点(即小正边形的顶点)在图像区域内的位置，确定其坐标。不同图像区域上的网格具有相同的结构和密度。附图12为图像区域上构造出的地形网格的示例。
7、确定地形网格点(简称网格点)的高程值：对于任一网格点p，首先确定它所属于的图像区域子块，记为B₀，并令图像区域子块B₀的相邻子块为B_i(i＝1，2，…，8)。从这9个子块B₀，B₁，B₂，…，B₈中寻找控制域覆盖网格点p的所有控制点，并把它们记为v_i(i＝1，2，…，k)，然后根据控制点v_i的高程值h_i、v_i的权w_i(p)、以及用户输入的全局控制参数η等计算p的高程值H(p)，如
H(p)=η1kΣi=1kwi(p)hi.]]>
计算出所有网格点的高程值后，即得到相应图像区域的地形块DEM。
8、建立大规模地形的DEM模型：每块图像区域上定义的DEM称为地形块DEM，上述得到的每个地形块DEM中的网格点的坐标都是局部定义的。本步骤分别为每个地形块DEM各设定一个偏移坐标，用以确定该地形块DEM在全局地形坐标系中的位置。各偏移坐标分别为相应图像区域最左下网格点的全局坐标。如附图13所示，对于整块图像区域的地形块DEM(i，j)的偏移坐标为(i×d，j×d)，其中d为正方形地形块DEM的边长。将地形块DEM的偏移坐标与块内网格点的局部坐标相加，便可得到该网格点的全局坐标。因此，将每个地形块DEM的块编号(i，j)、边长d、网格点的局部坐标等存储在磁盘上，并让相邻地形块DEM的边界网格点重合，便可构成大规模DEM地形模型数据集。
本发明方法相对于现有技术，具有如下优点：
1.提供多种参数对地形的生成过程进行控制，既能控制地形中局部细节的分辨率，又能控制地形块DEM的大小以及全局地形的规模。
2.生成大规模地形DEM模型数据的过程具有局部性，即：1)每个地形块DEM的产生过程仅涉及局部少量的标准化图像，有利于通过并行方法实现；2)单个地形块的产生过程与将要产生的全局地形的总规模无关。
3.可以生成任意规模的数字化虚拟地形DEM数据；能够通过在当前全局地形的边界处添加新地形块DEM的方式，对现存大型地形模型进行进一步扩充。
4.运用控制函数自然地实现了相邻地形块之间的无缝拼接，使相邻地形块之间的过渡更加自然，没有人工加工的痕迹。
5.对作为数据源的图像的选择要求不高，可以通过不同的组合方式由相同的图象集产生外观不同的多态地形。
附图说明
图1是本发明方法的流程图及相关数据的变换情况；
图2是经标准化处理后得到的正方形Lena图标准图像；
图3是对Lena标准图像进行半色调处理后得到的平面点集(脸部局部放大)；
图4是取K＝8时的平面点集四叉树层次剖分示意图；
图5是对图3中四叉树的叶结点进行点合并操作后的结果；
图6是由图3中的平面点集得到的1级简化点集；
图7是对1级简化点集再次进行四叉树层次剖分得到的结果(K＝8)；
图8是对图6中四叉树的叶结点进行点合并操作后的结果；
图9是由图5中的1级简化点集经再次简化得到的2级简化点集；
图10是控制点v的控制域，是以v为圆心、以R为半径的圆；
图11是一个控制函数的示例；
图12是在图像区域上构造地形网格的示意图；
图13确定图像区域的相邻关系，其中每个小四边形形表示一块图像区域；
图14是将每块图像区域一分为四，得到4块大小相等、形状相同的子块；
图15是由A，B，C，D四个图像区域构造的相邻区域；
图16是子块A_RD及其邻接的8个子块的示意图；
图17是相邻四块图像区域决定位于其中心的4个子块的相邻子块；
图18是图像区域的层次分割示意图及其四叉树表示。
具体实施方式
下面结合附图和实施例，对本发明做详细说明。
如图1所示，基于图像集的大规模数字化虚拟地形自动生成方法，其具体实施步骤如下：
根据用户输入的系列图像∑＝{I_i，j}以及地形制备控制参数(包括标准化图像的大小size、点集简化的层级level与点集分割控制参数K、地形网格的密度density、控制域控制域参数radius或num、地形起伏度控制参数η等)，构造大规模数字化虚拟地形。
1.对系列图像进行标准化处理：本实施例采用双线性插值法对大小为length×width像素的图像I进行缩放处理，将其转化为大小为size×size像素的正方形标准图像：首先以图像的左下角像素为原点建立直角坐标系，设标准化后的图像中任一像素为p(x，y)(x∈[0，size-1]，y∈[0，size-1])，通过变换得到其在原图像中的对应位置为q(x′，y′)，其中x′＝x·length/size，y′＝y·width/size，即有q(i+u，j+v)(i，j为正整数，u，v∈[0，1))。进而取p的颜色值为原图像中离像素点q′(i，j)最相邻的4个像素点a(i，j)、b(i+1，j)、c(i，j+1)、d(i+1，j+1)的颜色值C_a，C_b，C_c，C_d的混合，混合结果按照4个像素点与q′之间距离进行加权得到，即p的颜色值为
C_p＝(1-u)(1-v)C_a+u(1-v)C_b+(1-u)vC_c+uvC_d
注意到，对于原图像中的边界像素点，尽管与其相邻的像素点数不是4个，但上式同样适用，仅需将不在原图像中的像素点的颜色值设置为0即可。附图2示出了著名的Lena图经标准化处理后的结果。
2.对标准化图像进行半色调处理：本实施例采用Ostromoukhov提出的误差扩散算法来实现标准化图像的半色调处理，具体做法请参见：Ostromoukhov V.ASimple and Efficient Error-Diffusion Algorithm.In：Proc.ofSIGGRAPH’2001，2001，pp.567-572)。由此，得到所有图像的平面点集。附图3是对Lena图进行半色调处理后，所得到的平面点集的脸部局部区域的呈现情况。
3.简化平面点集：本实施例采用平面区域四叉树分割法实现平面点集的简化。如附图4所示，首先将图像区域分成大小和形状完全相同的4块；如果每块中包含的点的数目大于K(附图中取K＝8)，则将该块再次分为大小和形状完全相同、但尺寸更小的4块；对新产生的块迭代执行上述分割操作，直到每个新块中的点数都不大于K为止。整个平面点集的分割过程用一棵四叉树数据结构来表示(参见附图18)，最终不能再进行分割的块都位于四叉树的叶结点上。对于每个叶结点，如果其内包含的点数为N(N＞1，N≤K)，则在相应分割块内产生一个新点v来替代该块内的所有点v_i (i＝1，2，…，N)。其中v的坐标通过下式计算：
v=1NΣi=1Nvi]]>
这样就使四叉树中的所有叶结点块中只包含0或者1个点(参见附图5)。本实施例将此时四叉树中的所有点构成的点集称为1级简化点集(参见附图6)。
进一步地，删除当前四叉树中的所有叶结点，则四叉树中的部分内部结点变成叶结点。再次使用与前述方法完全类似的方法，将各叶结点块中点数多于1个的点合并为一个新点(参见附图7和附图8)，得到原平面点集的2级简化点集(参见附图9)。如此反复，直到点集简化级数达到用户要求的level级为止。最终得到的点集称为简化点集，简化点集中的点称为控制点。
附图4至附图9给出了平面点集四叉树分层简化过程的示例。
采用随机方法，为简化点集中的每个控制点v_i设置一个高度值h_i∈[-1，1]。
4.确定相邻子块并分解标准化图像：假定用户输入系列图像集∑＝{I_i，j}时，已经指定了各I_i，j之间的上下和左右相邻关系，例如：图像I_i，j与I_i±1，j、I_i，j±1以及I_i±1，j±1的相邻关系如附图13所示；若用户未指定上述上述相邻关系，则由计算机随机确定相关图像之间的相邻关系。附图13示出了一个平面四边形区域的图像区域铺嵌方式，其中每个I_i，j表示一块图像区域，而(i，j)表示I_i，j在铺嵌区域内的位置坐标。
将每块图像区域I_i，j分解为如附图14所示的4个均等子块，并确定每个子块的相邻区域，I_i，j的各子块的相邻子块在I_i，j、I_i±1，j、I_i，j±1以及I_i±1，j±1中。如附图15所示，有A、B、C、D四个相邻的图像区域，其中A中子块A_RD的8个相邻子块分别为A_LD，A_LU，A_RU，B_LU，B_LD，D_LU，C_RU和C_LU(附图16)。对于其它子块，可以类似地确定它们的相邻子块。在附图17中，相邻四块图像区域完全确定了位于其中心的4个子块的相邻子块，本实施例通过这种模式快速确定所有子块间的相邻关系。
5.确定控制点的控制域和控制函数：有两种方法来确定简化点集内每个点v_i的控制域：
1)如果输入的控制域参数为num，则通过寻找距离控制点v_i最近的其它num个控制点v_jⁱ(j＝1，2，…，num)来自适应地决定v_i的控制域半径。搜索各个v_jⁱ的过程需要遍历v_i所在图像区域及其相邻图像区域内的所有控制点。这样，便可取v_i的控制域为以v_i为中心、且恰好覆盖上述num个点的面积最小的圆形区域。即，控制域半径为：
Ri=maxj=1Num(||vi-vij||)]]>
2)如果输入的控制域参数为radius，则以数值radius作为简化点集内各点的控制域半径R，使所有控制点都具有大小相同的圆形控制域，即R_i＝radius。
本实施例采用有限支撑径向函数作为控制点的控制函数，即
wi(v)=[(Ri-||v-vi||)+Ri||v-vi||]2,]]>其中(a)+=a,a>00,otherwise]]>
它决定了点v_i的控制域中任意一点v相对于v_i的权重w_i(v)。
6.构造地形栅格：将图像区域的四条边均分为density等份，并称每份为段元；连接垂直和水平边上相应的分割点，可以得到如附图12所示的栅格网。设p(i，j)为栅格网中的任意网格点，则平面坐标(i，j)为p在X和Y轴方向上相对于最左下网格点所偏离的段元数目。如此，得到所有网格点的坐标。
7.确定地形网格点高程值：对于网格点p(i，j)，它必属于图像区域Ω的4个子块之一：Ω_LU，Ω_LD，Ω_RU或Ω_RD。不失一般性，令该子块为Ω₀，通过第4步所确立的相邻关系，可得它的8个相邻子块为Ω_i(i＝1，2，…，8)。遍历Ω_i(i＝0，1，…，8)中所有控制点，寻找那些控制域覆盖p的控制点，并把它们记为v_i(i＝1，2，…，N)，则可按如下公式计算网格点p的高程值H(p)：
H(P)=η1NΣi=1Nwi(p)hi]]>
其中，h_i和w_i(·)分别为控制点v_i的高程值和控制函数；η为用户输入的地形起伏控制参数，用于控制地形的起伏高度；
8.建立大规模地形DEM：由于每块地形的边长为size，所以图像区域I_i，j所对应的地形块DEM的偏移坐标为
             (i×size，j×size)
将地形块DEM的偏移坐标与块内各网格点的局部坐标相加，便可得到各网格点的全局坐标。让相邻地形块DEM的边界网格点重合，便可构成一个整体的大规模DEM数据集。