地形高程的高保真数字建模方法

专利名称：：地形高程的高保真数字建模方法
技术领域：
：本发明涉及一种地形数字高程模型的形成技术，具体说是一种基于局部面匹配地形数据采集模式和无缝嵌入矢量化技术的地形高程的高保真数字建模方法。
背景技术：
：数字高程模型产品(DigitalElevationModels,DEM)是囯家空间基础地理数据的核心数据产品，是GIS进行地形模拟分析和地学分析的核心系统，其应用遍布测绘、交通、军事、水利、农业、环境、资源管理、规划与旅游等众多领域，如在测绘部门主要用于三维地形建模、可视域分析、工程土方估算、遥感影像的几何矫正等；在地学界主要用于水文分析、蓄洪计算、淹没分析等三维地学分析；在国土资源部门主要用于土地利用类型详察等。构建国家级数字高程模型产品的两大关键问题是如何高效采集精确的地形特征数据和设计能适应复杂多变实际地形的数字高程模型的数据组织方式。现有的数字地形数据采集主要是釆用离散散点釆集和单点线串釆集两种模式。离散散点采集就是直接进行离散单点釆样，不釆集釆样点与釆样点之间的连接关系信息，如基于地形图采集离散高程点、野外釆集离散高程点、激光雷达扫描釆集等。单点线串采集就是以连续坐标串的方式采集和组织数据，如等高线的采集、数字摄影测量中地性线的采集，等等。上述釆样的共同特点是只釆集采样点的空间位置信息，不釆集釆样点作用区域的地形形态属性信息。实际上，地形表面的每一个采样点，它不仅是一个空间定位数据，更重要的是其形态属性信息，如对于地形特征线而言，在切线方向就具有是否是光滑连续曲线还是折线，在纵剖面方向是否可以进行光滑处理，光滑程度是多少等一系列属性。在理论上，每个釆样点局部范围的地形形态可以通过联合其周边采样点进行重构，然而在实际应用中，由于地形的复杂程度远远超出了现有数学与计算机的模拟能力，因此当釆样点局部形态属性信息缺失时，尤其是特征地形的形态属性信息缺失时，实际中很难高保真重构其地形曲面。事实上这也是造成现有DEM地形表达存在严重区域性失真现象的主要原因之一。现有DEM采用的数据结构格网DEM由于数据结构简单，易于计算机处理，且与遥感影像具有天然的相合性，因此获得最为广泛的应用，从目前的发展趋势看，DEM已经成为规则格网DEM的代称。目前，世界上许多发达国家如美国、英国、日本等不仅建立了覆盖本国的多种尺度的DEM,同时也开始推进高精度、高分辨率DEM，以及全球尺度SRTMDEM的建设，DEM成为国家经济发展、全球战略实施的核心支持系统。我国也已经完成国家范围1:1000000、1:500000、1:250000、1:50000和部分地区1:10000比例尺DEM的建设，在国民经济建设和科学研究中发挥着曰益重要的作用。大量应用实例表明，格网DEM虽能较好地反映山地和丘陵地区的地形自然起伏特征，但对自然突变地形、人工修整地形的模拟和存在严重的区域性失真问题。格网DEM数据结构的局限性和DEM构建技术的不完善是造成DEM地形模拟存在区域性性失真的根源。对于格网DEM数据结构的局限性问题，Ebner(1988年)提出在一般地区釆用规则格网DEM结构，沿地形特征处附加TIN数据结构的Grid-TIN混合结构。随后，许多学者对矢栅混合结构或矢栅一体化结构展开广泛研究(龚健雅，1992,1997;李德仁，1997;常燕卿，1998;杨树强等，1998;Amor,2004，2007)，出现了一些具有重要影响的研究成果，遗憾的是Grid-TIN结构的数据组织、维护和应用中的一些关键技术至今没有得到有效解决，还难以应用在大规模地形建模。龚健雅等1993年以格网数据模型为基础，提出多级格网数据结构，该模型利用在粗格网上进一步进行格网细分技术，实现对局部精细地形的描述，通过线性四叉树实现对多级格网数据的高效索引。多级格网数据结构对于二维平面信息的表达具有显著优势，但在形态高保真DEM中如何准确解算各个格网点高程值，正确维护精细地形的形态特征成为制约其推广应用的瓶颈。德国斯图加特大学的30^++程序将断裂线引入了格网DEM,实现了格网DEM对地形突变部位形态的精细描述，但其使用的数据结构、数据组织方式、数字地形分析模型及应用情况均鲜有资料可参阅，因此该模型还无法直接应用于实际DEM建设。同时，基于雕刻算法、距离函数法、表面生长法[朱庆，2006]、映射算法的三维空间数据模型也受到广泛的关注与研究，取得了一些重要研究成果。由于部分关键技术的制约，上述新型DEM数据结构都难以成为标准化的DEM数据结构，展开工程化建设。以上所述的现有DEM结构的缺陷，在精细化地形分析和工程应用中带来很大困难，在诸如平原地区洪水灾害模拟与防治、黄土高原水土保持与防治等工程应用上甚至会造成相关部门的决策失误，从而给国家和人民造成巨大的经济损失。
发明内容本发明针对以上缺陷，设计了一套新的高保真数字高程模型构建的工业技术方案。该技术流程涵盖了具有局部地形自适应匹配功能的推扫式地形数据采集技术和基于矢量特征数据无缝嵌入技术的格网DEM数据结构修正技术两大关键技术环节，该技术方案能够解决现有格网DEM地形模拟的区域性失真问题。本发明的基本思路是在不增加数据釆集工作量与工作强度的前提下，改变传统地形数据釆集的点釆集模式为局部面匹配釆集模式，完成釆样点空间位置信息和局地曲线形态属性信息的一体化采集，提高了地形信息釆集的准确性和完整性；并且，在保证原有格网DEM数据结构完整性和简洁高效性同时，在原有格网DEM数据结构中无缝嵌入矢量化的特征地形数据，从而实现地形的高效高保真数字模拟。具有广阔的工业应用前景。本发明地形高程的高保真数字建模方法技术方案如下第一步、推扫式数据釆集考虑到上述问题，本发明设计了针对复杂特征地形(见表l-2)的推扫式数据釆集技术。推扫式数据釆集的基本原理就是用一个具有一定空间形态的面状釆集器代替传统点状釆集器。在数据釆集过程中，推扫式数据釆集的面状采集器与采样点局部区域的地形精确匹配，在釆集釆样点空间位置信息的同时，自动识别和记录釆样点局部地形的形态属性信息。表l特征地形(PLA)内容与分类<formula>complexformulaseeoriginaldocumentpage8</formula>表2特征地形(PLA)空间形态属性表<table>complextableseeoriginaldocumentpage9</column></row><table>推扫式数据釆集的具体实现过程为(l.a)数据准备准备三维立体影像模型；(l.b)初始化釆集器(初始化采集器为现有技术)将初始化采集器处于面状釆集状态；该过程由推扫式数据釆集器的初始化函数Sample.InitializUon面状采集状态IntAreaSample完成,其中AreaSample是面状釆集的标识符；附表3介绍了推扫式数据釆集器的形态控制参数。(l.c)起始点釆集釆集起始点，此时采集器启动影像匹配线程Sample.TerrainMatching,等待用户釆集下一节点；(l.d)节点釆集依次采集其他采样点；当采集第一个节点时，Sample.TerrainMatching(影像匹配线程)依据起始点和相邻节点空间位置确定釆集器主方向，完成起始釆样点局部地形匹配，记录起始点形态属性信息；在随后的其他节点的釆集过程中，Sample.TerrainMatching自动依据当前釆样点和其最邻近的已经釆集的两节点的空间几何关系，确定前一节点的采集器主方向，完成其局部地形的匹配，获取其空间位置及形态属性信息；(l.e)终点釆集采集终点，采集器启动Sample.EndPoint—Matching(终点匹配线程)过程，依据终点和其最邻近两节点的空间关系推算终点位置釆集器主方向，完成局地地形的匹配；(l.f)非特征地形区域地形数据的釆集采集器调用Sample.Initializtion点状釆集状态IntPointSample过程，其中PointSample是点状釆集的标识符，设置采集器状态为点状釆集方式；然后按照现有技术的常规数据釆样模式采集所需的辅助地形数据；Sample.TerrainMatching和Sample.EndPointMatching的工作基本原理是依据三维立体模型中像元值与高程的相关性、相邻像元值变化程度完成地形匹配。如果三维影像质量不够理想，釆集器匹配精度不能满足要求时，就需要用户进行人机交互微调操作釆集器形态。表3是控制推扫式数据采集器的形态的参数说明。表4介绍了特征地形数据(PLA)空间位置数据的统一记录格式，表5介绍了标准化记录格式的特征地形(PLA)类对象结构，表6介绍了位置明显方式的PLA数据文件格式。表3推扫式数据釆集器的形态控制参数<table>complextableseeoriginaldocumentpage11</column></row><table>表4PLA空间位置数据的统一记录格式<table>complextableseeoriginaldocumentpage12</column></row><table>注为了标明坡面的坡度与坡向，倾斜面状PLA面内必须采集4个点，4个点应尽可能均匀布设在面域内，不可同处于一条直线。三个点用于计算坡度与坡向，一个点用于检查。表5标准化记录格式的PLA类对象结构〃PLA的数据结构StructsPLAIntPointsCounter.，IntCode^IntAttribution;DoubleRubdowir's3DPointPonits[];Doublez2[];〃PLA节点平面位置坐标数据结构Structs3DPointIntTag;〃节点属性，l表示原始坐标点，0为格网分割点。Doublex，y,z;〃位置坐标(x,y，z}//PLA节点个数，取l表示是特征点〃PLA的地形属性代码，备用属性，F-DEM构建并不需要〃PLA空间属性〃PLA磨光信息，强度或宽度〃PLA节点三维空间坐标〃面域的姿态信息表6位置明显方式的PLA数据文件格式〃PLA的数据结构信息ncols****nrows**"xllcorner****yllcorner****cellsizeNODATA—value""〃PLA数据列号行号PLA个数(n)PLA1的C0DE值形态属性编码节点数坐标信息PLA2的C0DE值形态属性编码节点数坐标信息PLAn的CODE值形态属性编码节点数坐标信息列号行号PLA的个数(n)PLA1的C0DE值形态属性编码节点数坐标信息PLA2的C0DE值形态属性编码节点数坐标信息PLAn的CODE值形态属性编码节点数坐标信息第二步、将釆集的矢量地形特征数据进行格网化分割就是用已有的格网DEM的格网线分割原始矢量地形特征数据，然后由DEM格网的四个格网点和格网内的矢量地形特征数据一起形成完整的数据组织结构；下文中的矢量地形特征数据均用字母PLA缩写代替。第二步所述的具体实现步骤如下(2.a)读取当前PLA的起始点作为漂浮点，计算并记录漂浮点所在的DEM格网的行列号，公式为Co/=,((浙-JTo)/Cs)肠=餅((取-yo)/Cs)(2.b)读取下一节点作为锚点，计算并记录锚点所在的DEM格网的行列号，计算公式同式(1);(2.c)判断漂浮点和锚点是否在同一个格网；如果是，且锚点是PLA的终点，跳转到步骤(2.f),否则，跳转到步骤(2.b);否，则进入下一步；(2.d)连接漂浮点到锚点，形成有向线PiPi+1(如图1所示)；计算水平格网线到有向线PiPw的角度6(水平格网线逆时针旋转到PiPw的角度，图l中的6角);计算并记录PiPw与DEM格网线的交点P(X，Y，Z)，计算公式如表7:表7格网化分割算法使用的主要公式<formula>complexformulaseeoriginaldocumentpage15</formula>'p指pi点到p点的空间距离，丄萍+,指pi点到pi+1点的空间距离(2.e)以步骤2.d中计算的P点作为新的漂浮点，跳转到(2.c)步骤，进行新的分割点计算；U.f)是否存在没有处理的PLA，如果是，顺序选择一个没有处理的PLA作为当前PLA，跳转到步骤(2.a);否，结東PLA格网化分割处理。上述处理过程全部集成在FDEM.PLAGridClip模块中，由计算机自动完成。第三步、将PLA进行类DEM格网化数据组织，经过格网化分割处理的PLA数据，按照栅格数据结构进行组织管理；图6说明了PLA与DEM—体化组织的基本流程；第三步所述的具体实现步骤如下(3.a)PLA节点与DEM格网点重合时，修正DEM格网点高程值与PLA节点高程值相同，仅在当前格网局地对象中记录PLA节点信息；(3.b)PLA分割点与DEM格网点重合时，如果二者高程相差不超过当前DEM规定的高程精度限差，直接在当前格网局地对象中记录PLA分割点信息；如果超过DEM规定的高程精度限差，修正DEM格网点高程值，然后在当前格网局地对象中记录PLA分割点信息；(3.c)PLA节点或分割点与DEM格网线重合时，与格网线共线的格网局地对象均记录该PLA节点或分割点信息；上述处理过程全部集成在FDEM.PLAGrid组织模块中，由计算机自动完成。第四步、生成新的数字高程模型(F-DEM)文件经过矢量特征数据的无缝嵌入，需要对修正的数字高程模型以文件形式输出，形成可以直接使用的测绘产品；其数据模型的设计既考虑了后续基于数字高程模型的地学分析方便，也着眼于和国家现有的格网DEM的兼容。本发明设计并提出一种特征嵌入式数字高程模型数据结构，该数据结构由3个子文件构成，包括格网文件、特征标识文件和特征文件；格网文件为常规的规则格网数据结构，与现有的国家DEM数据结构完全一致；特征标识文件用于矢量栅格一体化组织的空间数据索引，其结构也与现有国家DEM数据结构一致，但只记录有特征存在的格网；特征文件记录包含于一个格网单元内部的地形特征信息，包括特征点、特征线；相比传统的DEM格网化数据组织，特征嵌入式数字高程模型数据结构将格网化数据组织和矢量化特征数据组织通过特征标识有效融合，从而为后期工业化数据建库和基于海量数据的工程、地学应用奠定了基础。表8为该数据模型的详细文件记录格式，图7为该数据结构的示意图。表8新型数字高程模型的详细文件记录格式_(一)、格网DEM文件同现有格网DEM数据记录。(二)、特征标识文件(1)数据头数据头部分与DEM相同，但总行数和列数减1;(2)数据体如果该格网单元存在特征线，则该个网单元(i,j)的记录值为ID=(i-1)*(col-l)+j;反之以-9999标识，说明该个网单元无任何特征存在。(3)特征文件.-特征文件记录一个格网单元内部的特征分布信息，特征线与格网边界交点及其高程通过线性内插方式确定。特征文件按以下结构组织格网ID,特征代码1XI，Yl，ZlX2，Y2,Z2Xn，Yn,Zn特征代码2XI，Yl,ZlX2，Y2,Z2Xn，Yn，Zn注此处不需进行索引构建，格网本身就是最好的索引。第二~四步搡作实现了矢量特征数据的无缝嵌入技术。完成了对现有数字高程模型的高保真修正。本发明的技术特点及有益效果(一)本发明提出的特征嵌入式DEM数据结构，通过把特征地形信息通过格网化分割处理，无缝嵌入到格网DEM内部，采用"栅格化组织，矢量化模拟"技术，在维护格网DEM数据结构完整性和高效性的同时，修正了格网DEM对特征地形难以精细描述的局限性，实现了"全局高效，特征精细"的形态高保真DEM数据组织与管理，具有广阔的工业化应用前景。(二)本发明提出的推扫式数据釆集技术，在采集和编辑特征地形空间位置信息的同时，记录了特征地形在切线方向和纵剖面方向的形态属性信息，给后期DEM高程内插与地形表面的高保真重构与再现提供了强有力的数据支持，同时也降低了现有特征地形信息采集的工作总量，提升了数据采集质量，具有很好的实际应用价值。(三)本发明提出的F-DEM产品由三个文件构成，即格网DEM文件、PLA特征标识文件和地形特征(PLA)文件。三个文件可以独立使用、也可以集成使用。如果格网DEM文件、PLA标识文件、PLA文件同时存在，则执行F-DEM的高保真地学分析和精细化工程计算；如果只有格网DEM文件存在，则执行与常规格网DEM的地学分析与工程应用；如果只有PLA标示文件和PLA文件，或者只有PLA文件，则可以对PLA数据进行不规则三角网(TIN)的构建搡作，然后基于TIN模型进行工程和地学分析应用。图l:格网化分割算法使用的主要公式配图说明；图2:采集器在不同形态地形区域的俯视图(其中a为标准形状、b为传统形状、c为右侧为弯道时形状、d为左恻为弯道时形状、e为左侧为平直面时形状、f为右侧为平直面时形状)；图3:基于地形特征无缝嵌入方法的格网数字高程模型的修正流程其中A为一体化集成流程图(图中空心圆是PLA与DEM格网的交点、带x的圆是PLA特征点)B为PLA格网化分割示意图(图中空心圆是PLA与DEM格网的交点，实心圆是PLA的采样点。PLA格网化分割处理就是把PLA按DEM格网进行分块，每一个网内PLA按一个局地格网对象类进行数据组织管理。PLA格网分割处理的基本原则格网内PLA与相邻格网PLA的数据无关；不损失PLA的采样精度)；图4:本发明方法的流程图5:推扫式数据采集器在不同采样模式中的纵剖面图(其中a为三点法采样、b为五点法采样、c为两点法采样、d为单点法釆样)；图6:PLA与DEM—体化组织的基本流程；图7:修正后的新型数字高程模型文件数据结构图(其中a为格网文件、b为特征标识文件、c为特征文件)；图8:平原微丘地区由推扫式数据采集器采集的矢量特征数据文件；图9:平原微丘地区格网数字高程模型产品；图10:经过本专利修正算法修正过的数字高程模型产品(三维线框透视图)；图ll:经过本专利修正算法修正过的数字高程模型产品(分层效果图)。具体实施例方式下面结合附图和具体实施例，对本发明作进一步详细说明。(以平原微丘地形为例)平原微丘地形的总体特征是大范围的平坦地形嵌入了许多比较精细的规则地形，如围堰、导水渠、提坝、池塘、稻田、人工平直坡面，等等。国家现有的1:10000比例尺DEM规定的采样间距(格网大小)是5米，这在平原微丘地形一方面存在大量的冗余数据，另一方面难以准确重构和再现围堰、导水渠、堤坝、池塘、稻田、人工平直坡面等平原微丘地形的典型特征，造成现有DEM在该地区难以展开有效的地学分析与工程应用。该地区F-DEM的构建过程为第一步、推扫式数据釆集(一)数据准备采用航空摄影测量获得影像数据，在JX4中建立三维立体模型。(二)初始化采集器启动本发明的推扫式地形数据采集系统，初始化采集器，设置控制采集器的快捷操作键，做好数据采集准备；(三)特征地形数据采集设置采集器处于面状数据釆集模式，釆集所有特征地形信息；(1)起始点采集(2)节点采集(3)终点采集(四)非特征地形区域地形数据的采集设置采集器处于点状数据采集模式，依据实际情况采集非特征地形区域的点、线数据。第二步、调用FDEM.PLAGridCUp过程，完成PLA数据的格网化分割处理。第三步、调用FDEM.PUGrid组织过程，完成PLA数据的类DEM格网化数据组织。第四步、调用FDEM.Create过程，生成新型数字高程模型一一F-DEM数据文件经过该修正方法修正后的数字高程模型产品的三维线框透视图和分层效果图如图IO、ll所示。可以看出，相比较格网DEM而言，诸如提坝、水渠、梯田等特征信息完全表现出来，且数据精度高，保证了地形信息的高质J基于这种数字高程模型的洪水淹没分析等水文分析的精度大幅度提高。权利要求1、一种地形高程的高保真数字建模方法，其步骤如下第一步、推扫式数据采集(1.a)数据准备准备三维立体影像模型；(1.b)初始化采集器将初始化采集器处于面状采集状态；(1.c)起始点采集采集起始点，此时采集器启动影像匹配线程Sample.TerrainMatching，等待用户采集下一节点；(1.d)节点采集依次采集其他采样点；当采集第一个节点时，Sample.TerrainMatching依据起始点和相邻节点空间位置确定采集器主方向，完成起始采样点局部地形匹配，记录起始点形态属性信息；在随后的其他节点的采集过程中，Sample.TerrainMatching自动依据当前采样点和其最邻近的已经采集的两节点的空间几何关系，确定前一节点的采集器主方向，完成其局部地形的匹配，获取其空间位置及形态属性信息；(1.e)终点采集采集终点，采集器启动Sample.EndPoint_Matching过程，依据终点和其最邻近两节点的空间关系推算终点位置采集器主方向，完成局地地形的匹配；(1.f)非特征地形区域地形数据的采集采集器调用Sample.Initializtion，设置采集器状态为点状采集状态；然后按照现有技术的常规数据采样模式采集所需的辅助地形数据；第二步、将采集的矢量地形特征数据进行格网化分割就是用已有的格网DEM的格网线分割原始矢量地形特征数据，然后由DEM格网的四个格网点和格网内的矢量地形特征数据一起形成完整的数据组织结构；其具体步骤如下(2.a)、读取当前PLA的起始点作为漂浮点，计算并记录漂浮点所在的DEM格网的行列号，公式为Col＝INT((Xpi-Xo)/Cs)(1)Row＝INT((Ypi-Yo)/Cs)(2.b)、读取下一节点作为锚点，计算并记录锚点所在的DEM格网的行列号，计算公式同式1；(2.c)、判断漂浮点和锚点是否在同一个格网；如果是，且锚点是PLA的终点，跳转到步骤2.f，否则以锚点作为漂浮点，跳转到步骤2.b；否，则进入下一步；(2.d)、连接漂浮点到锚点，形成一条有向线；计算水平格网线到有向线PiPi+1的角度θ；计算并记录PiPi+1与DEM格网线的交点P(X，Y，Z)；(2.e)、以步骤2.d中计算的P点作为新的漂浮点，跳转到步骤2.c，进行新的分割点计算；(2.f)是否存在没有处理的PLA，如果是，顺序选择一个没有处理的PLA作为当前PLA，跳转到步骤2.a；否，结束PLA格网化分割处理。第三步、将PLA进行类DEM格网化数据组织，经过格网化分割处理的PLA数据，按照栅格数据结构进行组织管理；其具体步骤如下(3.a)、PLA节点与DEM格网点重合时，修正DEM格网点高程值与PLA节点高程值相同，仅在当前格网局地对象中记录PLA节点信息；(3.b)、PLA分割点与DEM格网点重合时，如果二者高程相差不超过当前DEM规定的高程精度限差，直接在当前格网局地对象中记录PLA分割点信息；如果超过DEM规定的高程精度限差，修正DEM格网点高程值，然后在当前格网局地对象中记录PLA分割点信息；(3.c)、PLA节点或分割点与DEM格网线重合时，与格网线共线的格网局地对象均记录该PLA节点或分割点信息；第四步、生成新的数字高程模型文件，经过矢量特征数据的无缝嵌入；特征嵌入式数字高程模型数据结构由3个子文件构成，包括格网文件、特征标识文件和特征文件；格网文件为规则格网数据结构；特征标识文件用于矢量栅格一体化组织的空间数据索引，其结构与现有国家DEM数据结构一致，但只记录有特征存在的格网；特征文件记录包含于一个格网单元内部的地形特征信息，包括特征点、特征线。全文摘要本发明公开了一种地形高程的高保真数字建模方法，其步骤如下第一步推扫式数据采集；第二步将采集的矢量地形特征数据进行格网化分割；第三步将PLA进行类DEM格网化数据组织，经过格网化分割处理的PLA数据，按照栅格数据结构进行组织管理；第四步、生成新的数字高程模型文件，经过矢量特征数据的无缝嵌入；特征嵌入式数字高程模型数据结构由3个子文件构成，包括格网文件、特征标识文件和特征文件；格网文件为规则格网数据结构；特征标识文件用于矢量栅格一体化组织的空间数据索引，其结构与现有国家DEM数据结构一致，但只记录有特征存在的格网；特征文件记录包含于一个格网单元内部的地形特征信息，包括特征点、特征线。文档编号G06F17/30GK101344390SQ20081002197公开日2009年1月14日申请日期2008年8月19日优先权日2008年8月19日发明者刘学军,汤国安,春王,旸陶申请人:南京师范大学