主要技术特点

7.2.1　主要技术特点

(1)卫星遥感图像地形影响消除技术

地形消除卫星影像主要应用于可燃物解译、正立体卫星影像制作、可燃物载量定量评估等。地物场指的地物电磁波的反射比表达的随机过程，它是以空间位置和太阳高度角的函数L(x，y，z，θ)为变量的随机过程{Z(L)|L(x，y，z，θ)∈T}，为地形构成的几何面对电磁波反射的影响量，对于确定的对地观测数据，θ是常量，所以L(x，y，z)就是传感条件的影响函数，Z(L)就是地物场，一般，可以把Z(L)看成单波段遥感图像或多波段遥感图像。

地物场Z(L)的一个样本就是一个对地观测数据集，其本质上是一个矩阵集，所以，有时也常用A(i，j)、B(i，j)等来表示单波段图像，i=j=1，2，…，n，也简单的记为A、B等。由于Z(L)是与地面位置有关的函数，而人认识和识别地物时，需要的是地物的内在表达，即地物的反射率，是与地形即地的起伏无关的图像信息，但这往往是难以办到的。为什么遥感图像既能看到地物又能看到地形，就是因为Z(L)是以地理位置为变量的随机过程。

图7.1　地形复原前后

如图7.1，将地形影响消除，就是得到与地的起伏、与太阳高度角无关的图像，变化其本质就是把随机过程Z(L)变换为随机变量р，р与地形、太阳高度角等无关，其物理意义是地物反射率，或消除了地形影响的对地观测数据。

光谱复原、定标后使阴坡和阳坡的地物色彩表达具有一致性。阴坡、阳坡的地物没有太阳光照的影响，使同物异谱的变化得到了恢复。遥感图像基本成为了反射率指数图像，这为遥感图像的计算机自动识别，为定量遥感实现，奠定了技术基础，如图7.2、图7.3。

图7.2　TM图像地形复原前后

图7.3　Spot图像去除地形影响前后对比图

(2)相对辐射定标技术

假设р(t)是与时间有关的反射率图像，р(Path，Row)是与轨道有关的反射率图像。定标的本质就是定义灰度的起点与单位刻度，并使所有的与时间、轨道有关的反射率图像定义在该度量空间中，并且有可比性。为了与计算机的可视化表达相一致，使用最小亮度刻度为0，最大亮度为255的模式空间。

定标的目的是使不同时间、不同空间上的对地观测数据具有可比性。从而使得它们能用一致的模型进行定量遥感估计，使它们合成的彩色图像颜色一致，生成的对地观测矢量与数据库一致。

同轨道的卫星遥感图像有90%以上的面积重合，相临轨道的卫星图像有10%左右的面积重合。相对定标是利用不同时间、不同空间的卫星图像的重叠性，以一个图像为基础，利用重叠部分的像素建立回归模型，将另一幅图像的灰阶计算为该图像的灰阶的方法。

相对定标也能实现不同时间、不同数据源、不同分辨率的卫星遥感数据的定标。当然，数据的重叠性与数据的统计规律的建立，仍然是定标的基本依据。

①相对定标技术方法

假设需要将图像B定标于图像A，用A_i表示A中的第i波段，用表示B中的第i波段。调入精确纠正的两景卫星影像A，B。在面积重叠区域，采集样本点集{(YA_i1，XB_i1，XB_i2，……i=1，2，3…)}，i为波段标号。

绘出{(YA_i1，XB_i1)，(YA_i2，XB_i2)，(YA_i3，XB_i3)…}的散点图，依据其走向，确定回归模型的基本形式。

进行数据变换与回归分析，求出定标模型。

对第i波段，依像素进行光谱度量的定标计算。

②相对定标技术要求

相对定标模型的相关系数大于0.8。

用统计量F进行回归模型检验，保证线性关系模型在95%的可靠性下，达到显著水平。F= sr/［se/(n-2)］服从自由度为1与n-2的F分布，其中sr为回归平方和，se为剩余平方和。

如图7.4，相对定标解决了不同轨道、不同时间的遥感图像的可比性，为定量遥感、图斑识别奠定了基础。

图7.4　相对定标前后的TM卫星图像

(3)DEM提取的地性线控制纠正技术

图7.5　基于DEM提取地性线示意图(山脊线、沟谷线)

几何精纠正是卫星影像应用的基础。进行遥感影像的几何精纠正处理时，常常面临控制点难找、控制点不足、配准影像精度低的难题。为解决该问题，提出了一种山地区域基于DEM地性线的控制纠正新方法，该方法以数字地形模型DEM为无几何变形的控制基准来纠正卫星影像。文中阐述了提取沟谷、山脊、山峰和凹地区域的地性线的原理与算法，给出了山地区域基于地性线进行卫星图像几何精纠正实施步骤，进一步讨论了地性线提取、控制点采集存在的问题，以及解决问题的途径。实验结果表明，对于山地区域，地性线的空间数量数倍于水系、道路等常规地图层;地性线来源于DEM，其空间稳定性和可靠性更高，可以用于山地区域的卫星影像的严格控制纠正。用该方法进行几何纠正处理，几何误差能控制在一个像元的水平上。

山地区域的卫星影像可以明显的反应地形的立体特征、包含丰富的地性线信息;DEM数据中包含有该比例尺水平上的所有地物特征点线，如山脊、沟谷、山峰、凹地等;用地性线的分叉、转折点等特征结构点采集控制点，再利用一般的几何精纠正方法，即可完成卫星遥感影像的几何精纠正。地性线的提取，控制点文件采集，多项式纠正处理，误差估计是几何精纠正处理的主要过程。

山地地形复杂多样，山脊线、沟谷线表达了地形的起伏结构，山峰、凹地表达了地形的高低结构。利用栅格分析方法，设计水流流向模拟算法、流量累积算法、地物的高点识别算法，实现山脊线、沟谷线、山峰点、凹地点的提取。

①沟谷线的提取

基于地表物质向低处运动的水流模拟方法，利用单流向最大坡降法确定每个栅格点的流向;依据每个栅格点上的水流方向，可计算出每个栅格点的累计汇水量。所谓沟谷，就是地形上的汇水区域，确定阈值，当某点上的汇水量大于该阈值时，认为该点是沟谷线上的点。阙值依靠人机交互选取，这个阈值与DEM的比例尺及其栅格点大小有关，取值太大，可能会导致表达小地形的地性线的丢失，取值太小，很多低地会形成面状闭合线，而不是最低沟谷区域附近的中线，合适的阈值是处理与分析的一个关键步骤。基于上述原理的地性线的提取包括了如下几个步骤:

a.栅格点的流向的计算:如图7.6所示，分别用1，2，4，8，16，32，64，128代表栅格A流向正东、正东南、正南、正西南、正西、正西北、正北和正东北的8个方向，遍历DEM的每个栅格点，如果某个点周边的高程相等，则该点无流向，记为0，单流向最大坡降法就是依据公式(1)，计算每个栅格点的8流向坡降值，取最大坡降方向为该点的流向。

Drop=△z/distance×100

其中，Drop为相邻栅格的坡降值，△Z为相邻的栅格的高程差，distance为相邻栅格的距离。遍历每个栅格点可计算出DEM的坡降方向。

图7.6　栅格点A的8流向

b.栅格点的累计流量的计算:栅格点的累计流量指的是流入某个栅格点的栅格数的累计值。一个栅格点的流向只有一个方向，但只考虑这个栅格点的周边8像素的时候，流入这个栅格点的水流最多是8，计算一个栅格点的累计流量，还要对整个地图中流入该点的像素遍历，即沿流入方向，向海拔高的地方寻找流入点，直到无流入点时停止，之后，统计流入该点的所有像点，得到该点的累计流量。依据上面的原理，使用8叉树建模，整棵8叉树的根节点为需要计算累计流量的栅格点，第1级子节点有8个分叉，8个分叉分别对应对于该点的左上、左、左下、下、右下、右、右上、上共8个栅格点(见图7.7，栅格点的流入模型)，如果这8个子节点的某个流向是向着父节点，则该分叉存在，否则，该分叉截断;之后，对于第2级子节点，如前面方法建立自己的儿子节点，直到任何一级的任何子节点不存在儿子接点为止。之后，统计根节点所生长的所有8叉树上的节点，该数即为累计流量数。使用优化算法，可去除重复计算量。

图7.7　栅格点A的流入与8叉树的父节点与子节点的关系

c.阈值的选取:调入累计流量栅格数据，显示为单色地图，亮度高为汇水量大的区域，对沟谷地形做基本的目视判断，选定沟谷的典型区进行交互查询，确定沟谷的累计流量阈值。并对图像做第一次沟谷提取试算，算法为:遍历所有栅格点，当累计流量大于阈值，Mgorge=1，否则Mgorge=0。Mgorge即为栅格化的沟谷地图。对Mgorge进行视化显示，判断提取效果，如果效果差，对阈值进行相应调整，直到满意为止。

d.矢量化计算:使用一般的栅格到矢量的转换函数，将栅格化的流域区细化为线状的矢量线，即生成了沟谷特征线，如图7.5所示。之后，还要对矢量地图做适当的清洁处理，如去掉非常短小的断线等。

②山脊线的提取

将DEM进行反转运算，把海拔高的区域变为低地，海拔低的区域变为高地。选择恰当的M，利用变换式:

DEM-1=M-DEM

一般来说，M最好取原DEM的最高高程的2倍。将原DEM实现地形反转之后，对新的DEM-1，重复沟谷线提取的过程，即可得到山脊线，如图7.5。

③单体山峰的提取

一般来说，单体山峰是点状地性线。山峰的海拔高于周边地物的海拔，将原DEM进行平滑处理，能起到削平高地的作用。将处理后的新的数字高程模型与原DEM进行比较，有较大削平作用的点即为山峰点。使用一般重采样方法，对原DEM进行重采样，像素尺寸取原来的2倍，得到新的数字地形模型，用平滑滤波方法对其进行2到3次平滑。把原DEM与新的数字高程模型做差值运算，用(1)的③、④中的方法，计算出山峰点矢量数据。

④点状凹地的提取

用(2)做DEM的反转运算，用(3)的方法提取单点凹地。

在一般GIS系统的地图显示与编辑环境中，将经过系统级纠正的ETM卫星图像调入并显示为假彩色;调入沟谷线，用蓝色显示;调入山脊线，用红色表示;调入山峰地图，用红色的点显示;调入凹地地图，用红色点显示。手工交互采集纠正需要的沟谷控制点，即用同名地物点匹配方法，采集沟谷或凹地上的特征点以及其对应的卫星影像的同名地物点的坐标。将全部数据旋转180°后，采集山脊控制点，并生成控制点文件。

由于随机变形误差具有不均匀性、不一致性，卫星影像的每一小块区域周边都有自己的控制点，不能使用适用于整体区域的高次多项式或有理函数纠正，使用的是自适应三角形网络法纠正。例如，纠正前后水系与影像匹配对比效果见图7.8。

卫星图像几何纠正存在的误差是随机变量，依据误差理论，其服从正态分布，依据统计估计理论估计误差。误差估计的样本点按如下方法实现，建立一个公里格网地图层，公里格网的起点满足随机性要求。以格网点为单元进行随机抽样。对于抽中的格网点，寻找离该格点最近的可辨别的明显地物点，以地性线图层为基准，在计算机屏幕上的测量卫星图像还存在的几何误差，各个格点上的误差集合X₁，X₂，…，X_n构成了误差估计的独立同分布随机样本。用样本的均值ΔX估计总体误差。ΔX的估计式为:

这里，n为样本点数目。置信度为1-α的置信区间为:

［ΔX-St_α/2(n-1)/n^1/2，ΔX+St_α/2(n-1)/n^1/2］

这里，S²={(X₁-ΔX)²+(X₂-ΔX)²+…+(X_n-ΔX)²}/(n-1)，t_α/2(n-1)为t(n-1)分布的α双侧分位数。

图7.8　纠正前后水系与影像匹配对比图

(4)卫星遥感图像的正视觉立体影像制作技术

航天遥感摄影的机理和物理学基础决定了遥感成像的是反视觉的，即负立体图像，这增加了遥感图像在行业和区域应用推广的难度。将负立体卫星影像转换为正立体卫星影像，解决了基层林业单位人员不熟悉、不习惯原卫星图像固有的反视觉地形表达的问题。使用亚像元分解、图像增强技术，能使制作出的正立体影像具有更好的可视化效果，便于遥感图像在林业上的业务化推广应用。例如，图7.9中显示的基于TM影像制作的正立体影像图，将更加有利于集成到森林灾害监测、应急指挥GIS业务系统等业务化应用平台中。

图7.9　负视觉立体卫星图像与亚像元分解并增强后的正视觉立体卫星图像对比

(5)可燃物分布图斑的自动识别技术

该技术可用于森林可燃物类别识别、载量评估。在进行了光谱恢复及增强处理的多个遥感图像中，选取光谱复原较好的波段，组合成多波段彩色堆栈，从多光谱中提取出有光谱特征和空间特征区别的图斑，将识别出的图斑根据其值的不同进行矢量归并，转化成矢量多边形，栅格—矢量化时，图斑的光谱属性值即本征覆盖类型集成到属性库中。这时提取出的图像的矢量多边形，是表达矢量地图的最小信息载体单位。同时，对该矢量多边形中非常小的、破碎的多边形进行去除，该研究中去除了面积小于15亩的多边形，即小于11个像元的斑块将被去除，合并到相邻的大的图斑中，最后进行部分矢量边界进行手工修改与编辑。效果图如图7.10所示。

图7.10　自动“图斑识别”与多光谱分类矢量对比

图斑定量遥感因子的提取是图斑识别矢量化后进行语义解译时重要依据，特别是对于划分植被类型时，各定量遥感因子是植被类型的主要依据，部分可集成因子如下:

①图斑表面均匀性

纹理度量是描述图像的均匀性、一致性及其格局重复或相似性的度量。在目视解译时，往往根据纹理元素的度量作为依据。由于大脑思考时利用纹理的复杂性，所以，用计算机抽象定义图像纹理的指标也较为困难。本文在根据图斑的像素点的均方差，来描述图斑均匀性、光滑性的纹理度量指标。

Std2=sprt{1/n×∑［B(i，j)-mean］}

其中，n为图斑内所包含的像点总数，B(i，j)为图像在i行j列的灰度值，mean为图斑内所有像点的平均值。

②近红外反射率指标

由于绿色植被叶子的海绵组织对过0.8～1.3m的近红外光强烈的反射，形成植被反射光谱曲线上的最高峰，它直接对应于地面植被的叶面积指数，可用来估计植物的密度等。该指数因植物群落的类型、组成、健康状况的差异而呈现不同的值。

该指标主要对应于TM4通道。

③红光反射率指标

绿色植被由于内部细胞特性，在进行光合作用时吸收红光作为能量，故有较强的红光吸收，该指标对应TM3通道。能体现植物光合作用状况，用于估测蓄积量和生物量，植被健康状况，估测土壤裸露和土壤侵蚀等，进行植物或树种类型的智能化推理。

④归一化差值植被指数(NDVI)

归一化差值植被指数指数对于突出植被特征、提取植被类型或估算植被生物量。该指标的定义是根据植被在红光处的反射低谷与近红外处反射的波峰，即对红光与近红外波段进行组合比较运算。对于TM数据，它的运算为:NDVI=(TM4-TM3)/(TM4-TM3)。

⑤比值植被指数(RVI)

该指数与NDVI作用相似，与植被的生物量状况有密切关系，只不过它定义为近红外波段与红光波的比值运算。对于TM数据，它的运算为:RVI=TM4/TM3。

⑥汇水指数

汇水指数反映的是沿着水流方向，汇流到各像元的总的水分量。它是通过影响植被的水分平衡来影响植被的分布，定义为每个像元在其坡向上的汇水多少。

对以上求出的各定量遥感因子是基于各像元的，通过运算，将其集成在矢量多边形元的属性库中，如图7.11。可以通过设计利用矢量多边形中各像元定量遥感因子的平均值的算法，即利用GIS空间栅格选取功能，选取落入矢量多边形的像元，求其平均值作为矢量多边形相应的属性值。同时，也可以将生态区区划矢量图转换成grid栅格，通过与多边形矢量的空间运算，集成矢量多边形所属的生态区。

图7.11　经遥感定量估计因子集成后的遥感影像数据

(6)用户化二次解译技术

如图7.12，加硬拷贝界线(红)识别后的Spot图像、用户化解译为林地“小班”。用户基于指标和数据库属性定义用户需要的语义地类，通过数据库操作和多边形基元的合并操作，实现用户的二次解译。专业语义解释是提供给用户二次解译、使用对地观测矢量的一种技术方法，是在用户端实现的。该方法只提供基本原则和质量控制标准，用户可以充分依据从样本数据中发现的规律，来进行专业的语义描述。

①用户必须依据自己的专业应用需求，建立专业分类标准与分类系统。

②准备判读样地数据集(最好是GPS样点)。

③判读样地与对地观测矢量进行叠加分析，建立语义解译规则。将判读样地数据与对地观测矢量进行叠加分析，确定本征地物类到专业类的对应关系，发现某个类别所具有的个性的特征规律，充分利用定量遥感指标，建立描述个性规律的模型，建立利用属性数据库中的定量遥感指标参与语义推理的准则。

④用数据库语言实现语义解释和推理。

图7.12　加硬拷贝界线(红)识别后的Spot图像

(7)亚像元分解与可视化增强技术

利用亚像元分解与增强方法，增加了卫星图像的可读性、可视性，消除马赛克效应，如图7.13～7.16。卫星图像用可视化的方法，把直观、客观的地理信息传递给人，所以，图像的可视化增强就成了一项卫星图像可视化应用的基础工作。双线性内插值法、均值平滑和拉普拉斯变换形成的增强图像函数是数字图像显示处理的一般方法。双线性内插值法缩放后图像质量高，不会出现像素值不连续的情况，并由于其具有低通滤波器的性质，使高频分量受损，所以可能会使图像轮廓在一定程度上变得模糊。在获取的原始图像中一般带有一些噪声，为了消除这些噪声，均值平滑通过平滑滤波器，把突变点的灰度分散在相邻点中来达到平滑效果，操作起来也简单，但均值平滑也通常带来图像模糊，而且图像尺度N(领域)选取得越大，模糊越严重。即使应用了拉普拉斯变换形成的增强图像函数也不能弥补部分丢失的图像信息，处理后的图像质量并不十分理想，立方卷积内插法的效果类似。并且，这些内插法针对亚像元分解效果不佳。

亚像元分解专利技术采用的方法是直接对数字图像像元进行如下分解、回代与内插，可消除以上方法带来的弊端，得到可视化效果更好的图像。利用双线性内插算法对图像G1进行像元大小S/2的亚像元分解，即把一个像元分解为4个像元，得到图像G2;用G1的数据值替换G2中偶行偶列的数据值得到G3，即:①分别输入N×M数字阵列A数字图像G1和2N×2M数字阵列B数字图像G2;②判断G1的i、j大小，如果小于N、M，则将x=2×i、y=2×j赋予G2的x、y，并将阵列Aij的数值赋予阵列Bxy，合并所有的Bxy，得到图像元G3;③再利用双线性内插算法对G3进行像元大小S/3的亚像元分解，即将一个像元分解为9个像元，得到G4;④用G1的数据值替换G4中3×3区域中心位置的数据值得G5，即:a分别输入N×M数字阵列A数字图像G1和3N×3M数字阵列C数字图像G4，b判断G1的i、j大小，如果小于N、M，则将x=3×i、y=3×j赋予G4的x、y，并将阵列Aij的数值赋予阵列Cxy;⑤合并所有的Cxy，得到图像元G5;⑥对G5先进行均值平滑处理，再进行拉普拉斯增强处理;最后再用双线性内插算法对G5进行像元大小为所需值的处理。

以上步骤③、⑤像元分解的2规则和替换规则，可以重复地单独或组合使用，从而将像元分解为任意值。也可以将原图像的像元分解为任意大小，对像元分解与增强所形成的亚像元进行原图像数据值的部分替换，再应用通常的处理方法，得到所需要的较为理想的数字图像亚像元矩阵值，实现图像质量的增强。

图7.13　原始的30m分辨率的反视觉立体影像

图7.14　分解为5m分辨率的正视觉立体影像

图7.15　经过亚像元分解、增强处理后的安宁市陆地卫星图像

图7.16　地形图立体可视化处理结果