内蒙古草原牧场防护林遥感技术方法研究

内蒙古草原牧场防护林遥感技术方法研究_巡天察地梅安

内蒙古草原牧场防护林遥感技术方法研究

梅安新^[1]

摘　要：本文引用了多种航天遥感信息源，以光化学与计算机相结合的方法进行专题信息提取和自动分类。在此基础上，结合当地的特点提出了TM信息最佳波段组合，为解译提供了优质图像的处理流程，并为土地资源利用评价提供了自动评判程序，对内蒙古草原牧场区造林立地条件和土地合理使用进行了多目标预测。

关键词：防护林遥感，波段组合，自动评判，多目标预测

一、引　言

内蒙古草原牧场防护林区，共有35个重点县（旗），1个典型县，面积达22 000km²，相当于江苏省、浙江省两省面积的总和。要在三年时间内完成森林分布图、草地资源图、土地利用现状图、造林立地条件图和土地资源评价图等5种专题图和相应的面积数据及造林生态效益的分析研究，工作量大，依靠常规地面调查工作是无法完成的。在调查研究过程中，我们不仅以航天遥感影像作为主要信息源，而且还对遥感技术方法进行了探讨，提高了影像解译的效果和精度。

根据研究区遥感信息源的条件，我们尽可能地获取多种遥感资料，其中包括：Landsat TM CCT资料9景，TM影像31景，SPOT HRV CCT资源1景、中国国土卫星影像（覆盖库伦、奈曼、翁牛特及凉城等4个旗县），NOAA AVHRR CCT资料8景。此外，还有喀喇沁旗旺业甸林场600km²的彩红外航片，对航天遥感资料，我们选择了Landsat的TM作为主要信息源，进行了专题信息提取、计算机分类、动态分析，土地资源评价自动评判及区域多目标规划等研究。

二、遥感图像处理和信息提取

遥感图像处理的目的有二：其一是为目视解译提供可读性好的、主题信息突出的优质影像；其二是进行自动分类和动态、数量分析。对于大范围的区域研究来说，不可能每一景图像都进行数字处理，因此，我们以计算机和光化学相结合的方法为主，对于典型地段则进行数字处理。

1.TM图像处理

①TM影像最佳波段组合的方法与机理研究

目前使用最广的TM合成方案是Landsat TM_4（R）、TM_3（G）、TM_2（B）的合成，相当于MSS₄，从最大离差值、最大信息量及最小相关性等方面对波段最佳组合进行了研究，提出了MSS₅，MSS₇的合成，被称为标准段彩色合成。研究表明，这一组合并不是最佳的组合，远远没有充分发挥TM多波段丰富信息的潜力。我们从机理上、方法上研究结果表明，经过分段拉伸处理的TM_4（R）、TM_3（G）、TM_5（B）合成的图像解译效果最佳。最佳波段组合选择的原则有三，即：最大离散度、熵最大值和最小信息相关准则。

a.离散度最大准则

离散度表示各种地物的分离程度，离散度愈大，各类地物愈易区别。这里以欧氏距离作为分离性的测度。

设n维空间中二点Z1＝（Z₁₁，Z₁₂，…，Z_1n）和Z₂＝（Z_21，Z₂₂，…，Z_2n），则Z₁，Z₂的欧氏距离为：

d＝（Z₁－Z₂）′（Z₁－Z₂）

设各地物在影像上的密度矩阵为：

从（1）式中选出N×3的子阵：

其中a₁，a₂，a₃为不相同的任三个波段，则距离测试为：

其中：Z_i＝（X_ia₁ X_ia₂ X_ia₃）′；

　　　Z_j＝（x_ja₁ x_ja₂ x_ja₃）′。

取d值最大的（a₁，a₂，a₃）组合即为最优段组合。

b.熵最大准则

熵是信息量的一种测度，代表了信息量的相对大小，熵愈大，所含信息愈多，假定M个波段中选出的m个波段的联合的高斯分布的情形下，设矩阵（1）的协方差阵为：

经过推导，得出信息指标值：

I＝（X_a1，X_a2，…，X_am）＝｜C_m｜

C_m是（3）式的一个子协方差阵：

｜C_m｜是m个波段变量的子协方差阵C_m的先烈式值。当C_m的特征值已知时，信息指标可取另一形式：

其中λ₁，λ₂，…，λ_m是子协方差阵Cm的m个特征值，信息指标是子协方差阵的m个特征值的乘积。取I值最大的（a₁，a₂，a₃，…，a_m）组合即为最优组合。

c.信息相关最小准则

对于卫星影像来说，每景影像单波段的信息量可以表示为：

式中：H——影像信息量；

　　　P_i——第i级灰度出现的频率。

入选的每一个波段的信息量H值都大，其合成影像的总信息也不一定必须考虑各波段之间的信息相关系数R_ij，要得到R_ij只需在计算机中读入TM各波段数据，即可计算出各波段的信息相关矩阵，选其中三个相关最小的波段进行合成，即为最优波段组合，由于这方面研究，国内外都做了大量工作，不作详细讨论。在这里我们综合考虑信息量和各波段之间的相关系数，提出了一个波段选择函数：

式中：P（H，R）为选择函数；m为合成波段数目（在通常情况下取m＝3）；H_i为第i波段的信息量；R_ij为i波段与j波段影像之间的相关系数，直接利用P（H，R）值计算选择波段比较麻烦，因此，可以在使P（H，R）值保序的情况下，即不改变波段组合后P（H，R）值的大小时，作一近似计算，使此函数能很方便地进行计算。

对于H值来说，由于影像概率密度近似正态分布，则每一波段的信息量可写为：

其中：

式中：σ——方差；

　　　μ——均值。

由此可以导出：

L_nx是一种单值正交函数，故可省略L_n形式并舍去常数，得到一种相对信息量指标值：

式中：σ——概率密度曲线方差。

根据“3σ规则”，正态随机变量落在区间［μ－3σ，μ＋3σ］内几乎是肯定的。对于影像概率密度来说，相当于灰度出现的区间N_ε［n₁，n₂］，则

其中：n₁——影像灰度最小值；

　　　n₂——影像灰度最大值。

可见，每幅影像的信息量可以由影像灰度出现阶数近似地表示。

前面提到影像相关系数的计算是通过计算机图像处理系统来实现的。这里提供了一个可以利用胶片密度值测量，通过计算机计算的方法。其基本原理是利用同类地物测得波段间相关系数近似地代替CCT数据的波段信息相关。对于三个波段的合成，则

其中n₁，n₂的值可从胶片上测得，R_ij值可以从有关资料中查得。为了方便计算P（H，R）值，我们在计算机中编制了一个运行程序，只要输入n₁，n₂及R_ij值即可打印出P（H，P）最大的三个值及其相应的波段组合，其框图如图1所示。

我们对Landsat 122－30幅TM影像进行试验，除TM₆外，其余各波段在最优组合中出现的频率见图1。

图1　波段选择程序框图

其中TM₃、TM₄、TM₅出现的频率最高，与计算机上所做的最小相关系数计算的结果是一致的，同时，也与各波段信息的物理特征是一致的，分析表明，在TM的7个波段中，可以组合成亮度、绿度、湿度、热度、透射度等物理量。通常以亮度、绿度、湿度三个物理量最为常用，通过对上面三种方法所选择出来的最优波段组合TM₃、TM₄、TM₅，可以近似地表达上述三个物理量，TM₃可以近似地代表地物在可见光波段的亮度特征；TM₄可以近似代表绿度，表示植物生长状况；TM₅则可近似代表地面湿度状况。在湿润地区用TM₅代表湿度往往偏大，掩盖了湿度的差别，但对于半干旱和干旱地区来说，TM₅则可近似地代表地面的湿度状况。这三个波段的组合，信息量最丰富，相关性最小，而且能近似代表三个物理量。这个结论在实际处理119—28，120—29，120—30，121—29，121—30，121—31，122—29，122—30等幅影像中完全得到证实。以TM_4（R）、TM_3（G）、TM_5（B）经计算机三分之一拉伸和进行光学合成处理效果最好，各类地物在影像上有明显的差别。以120—30（1989年5月5日）为例，生长良好、茂盛的植物如大片林地等呈鲜红色，疏林地呈淡红色斑状，灌木林呈暗红色斑状；水浇地在暗绿色背景上有红色斑点；耕地中的林网在深绿色背景上呈红色网格，非常明显；旱地呈淡肉红色；沙地呈白色；盐碱地呈黄色；水体呈绿色；草地呈茶色，整幅影像色彩非常丰富，是TM_5（B）、TM_3（G）、TM_4（R）合成，以红蓝两色为主的图像所不能比拟的，这就为目视解译提供了可判程度很高的优质影像，大大提高了解译精度，开鲁县1∶20万解译图上的林地图斑与当地同比例森林分布图上的图斑符合率高达90%以上。

实践表明，通过机理研究和方法探讨，选择最优波段组合，为解译提供高质量影像，这是提高解译精度的基本保证。

②在几何配准的基础上，TM图像的空间卷积处理

许多研究已表明，TM影像光学放大时，1∶10万是最佳的比例尺，此时，影像的构像潜力得到充分的发挥，目视效果最好。但是，为了制图的需要，仍可进一步放大，在具体处理过程中，我们将赤峰市图像CCT数据，先分波段经两次卷积插值放大，然后再按所需波段进行合成，获得1∶5万彩色影像，甚至可放大到1∶2.5万，仍有较好的目视解译效果。经过多项式拟合几何校正的这类影像，可以满足1∶5万目视解译和专题制图的要求。

研究还表明，经过下列步骤处理的影像，可为目视解译提供最佳选择：

a.几何校正（多项式拟合或刚性平移处理，视所需精度而定）；

b.三次褶积；

c.分段线性拉伸；

d.Laplace变换；

e.最佳波段组合选择或特征空间变换处理。

③专题信息提取

主要根据本次遥感调查的要求和本地区的特点，我们利用TM图像光谱段多的优势，进行光学专题信息提取，其目的是突出目标的同时，压低背景噪声。

基本方法是通过对影像的分析，发现有些地物（主要是林地）在常规的合成中与农田上的植被很难区别，而在适当的时相，在某些波段中，两者有细微的差别，如赤峰幅1988年10月6日一景，奈曼幅1988年9月30日一景，在TM₄影像上有微小差别，我们在标准假彩色合成（TM_2＋B，TM_3＋G，TM_4＋R）的基础上，再以TM_4＋G曝光一次，叠加在标准假彩色影像上，使林地颜色成黄色（R＋G＋Y），与农田晚秋作物分开。此外，还把预先作出的TM₄/TM₃比值影像与其他波段合成，取得植被信息突出的影像。在处理过程中还通过控制曝光量等方法来增强我们所需目标信息，有利于林业、立地、土地利用等专题的解译。

④计算机分类方法探讨

目前公认的监督分类中精度最高的是最大似然分类法。但是，由于运算工作量很大，所耗机时多，使它的广泛应用受到影响。为了探索更经济、更快速的方法，我们对各种方法进行比较，选择了Kauth－Thomas变换，其变换矩阵如表1。

表1　Kauth－Thomas变换矩阵

通过对127—33幅伊金霍洛旗范围处理的结果表明Tasseled Cap变换图像与最大似然分类的结果基本符合，但其所耗机时仅需几分钟（相应图幅范围最大似然分类需几个小时）。因此，在某些特定的地点和条件下，可以用Tasseled Cap变换来代替最大似然分类，以节省大量机时。

2.SPOT图像分区分类处理

在研究区内，我们收集到一景1986年5月7日的SPOT卫星图像，先进行目视解译，分成若干自然区（如山区、平原等），分区聚类，划分类型。从而解决了山地和平原很难采取同一判别函数的问题，大大地减少了混分、错分现象，这一方法的实质是把人的知识与计算机图像处理结合起来，提高了分类精度，使原先只能划分6、7个植被单元增加到21个。

3.NOAA卫星的宏观动态分析

NOAA卫星由于其空间分辨率的限制（AVHRR的星下点分辨率为1.1km），不能满足调查制图要求，但由于其覆盖周期短，一天可获得两幅图像，对宏观动态监测特别有利。

在常规处理的NOAA图像上，除了科尔沁沙地的影像与周围有明显差别外，其沙地内部及外围的防护林信息都难以显示，我们这次把NOAA CCT数据进行了多种处理试验，得到两种效果较好的图像。

首先，我们求出“绿度”值（2/1），再与其他波段（例如Ch₃－3.55－3.93_μ）进行合成和伪彩色处理，或直方图均衡加伪彩色处理，使植被信息、地面湿度信息在背景上突出来，从而区分出植被较差的沙坨地和植被较好的甸子地，沙地中还可区分出流动沙地、半固定沙地和固定沙地等。特别是科尔沁沙地内部植被恢复得较好的吉尔戈郎、甘旗卡等地，大青沟自然保护区及沙地南缘的“三北”防护林立体工程和黄土丘陵台地上的人工林都得到明显的反应。说明“三北”防护林体系建设已取得明显的效果。

其次，还对不同时相的NOAA AVHRR图像作同样处理，例如将K－L变换的第一主分量进行分割，能较好地反映不同时相植被景观的结构。可以认为，通过K－L变换的新图像代表了植被景观的“平均”状态，不同时相的差异，即为其动态变化。

三、土地资源评价的自动评判与机助制图试验研究

土地资源评价长期以来都是以人的经验为主的评判来确定土地等级。但随着遥感技术的应用及信息系统学科的发展和计算机的普及，为土地资源评价的自动评判提供了有利的条件。我们在科尔沁左翼后旗的土地资源评价中，应用计算机发展了自动评判的软件。主要步骤和做法如下：

1.评价因素的选择

主要依据《“三北”防护林遥感综合调查技术规程》的规定，应用因素分析法选出主导限制因素。主导限制因素的选择和权重计算，是通过计算原始数据的内积矩阵的特征根与特征向量实现的，特征根的大小反映了主因素包含数量信息的多少，而特征向量则反映了主因素与变量之间的相关关系。每个主因素的意义不仅仅从初始载荷矩阵中反映，而要对因子载荷矩阵进行旋转，使其结构简化，以再现主因素和原始变量的内在联系（求解过程从略）。基于这一方法原理，在计算机上设计了程序并进行运算，得出水分（w）、细土层厚度（d）、盐碱化程度（s）、潜育层厚度（h）、土壤质地（m）、风积风蚀（v）作为主导因素。这一结果，对科尔沁沙地来说，与理论分析结论是一致的。

在科左后旗土地资源评价中我们应用了模糊综合评价方法。考虑到各质量因素的权重与单组质量因素指标的复合作用，也就是一个以f（u）到f（v）的模糊变换R，给定定权分配A∈f（v），则

B＝AOR

就是以A为权数分配的模糊综合评价。

评判中的因素集应该是同时倾向一个目标的，比如从各限制因素考虑土地的可耕程度（宜农）、风蚀风积等级越小越好，水分条件越好越可耕，都是因素集中的方向因素；对象的处理是各自进行的，相互间没有关系；评判函数简化后即是求最大、最小排序；权重的选择则可由主因素分析法中载荷因素值来确定。其综合评判步骤如下：

①给出评价单元集：x＝（x₁，x₂，…，x_n）；

②找出限制因素集：v＝（u′，u₂，…，u_m）；

③找出评判矩阵：h＝x×u→［0，1］

r_zj＝R（x_i，u_i）←［0，1］

其中，r_zj是评价单元x_i在限制因素u_i上的特性指标。显然，x_i的特征向量为R1_xi＝（ri₁，ri₂，…，ri_m），i∈［0，1］，故可视为V中的Fuzzy；

④确定评判函数f＝［0，1］^m→R，记

D＝f（Z₁，Z₂，…，Z_m）

其中D＝max（ri₁，ri₂，…，ri_m）（i≤m）。

这里的评判指标即为土地等级的限制强度级映射［0，1］^m。

最后，把求得结果在ATLAS GRAPHICS自动制图系统支持下，在计算机上分块绘制科左后旗土地资源评价图。并经MAPEDIT模块数字化获得的图形数据文件＊＊.BNA和DATAEDIT模块人机交互处理后的属性数据文件＊＊.DAT，采用AGAREA面积计算程序，求出科左后旗各等级土地的面积（见表2）。

表2　科左后旗各等级土地面积表

这一试验研究的结果表明，遥感解译的成果可以应用计算机进行土地等级自动评定的制图，代替大量人工作业，还可作为一子模块，加入GIS，扩充GIS的性能。

2.应用遥感的解译变量进行旗（县）最佳土地利用结构模型的研究

我们在进行“三北”防护林遥感调查的同时，还结合区域的开发和合理的利用土地资源，使“三北”防护林的建设与当地农牧林各业的生产发展协调一致，可使两者更紧密地结合起来。因此，在奈曼旗进行解译结果分析的同时，建立了多目标规划数学模型，为解决土地资源合理利用提供科学分析手段，也把“三北”防护林遥感调查成果的应用更加深化。

线性规划作为一种优化研究技术已较成熟，其应用也相当广泛，但也有其难以克服的缺点，主要是其目标单一以及资源的有限性和目标需求相矛盾，而在规划实践中，涉及的大都是多目标问题，特别是没有统一度量单位的多目标问题，线性规划更是无能为力。另外，在约束条件中出现互相矛盾的方程时，线性规划可行解域而出现无解的情况。我们在奈曼旗应用遥感作出土地利用图和土地资源评价图之后，把两者结合起来分析，并应用多目标规划的数学模型，用优先等级方法，提出最佳土地利用结构，即首先按各种目标的重要性排成顺序，在求解过程中，只有保证实现前一目标的前提下，才考虑下一目标，并且在求解后面目标时，不允许改变已得到的前面目标达成值。

在建立奈曼旗土地利用结构模型时，选择了14种土地利用方式作为规划决策变量。又为了便于计算畜牧业产值，把大畜、小畜、猪等也作为决策变量，共计17个决策变量，即：粮豆作物面积（x₁）、经济作物面积（x₂）、其他作物面积（x₃）、饲料地面积（x₄）、防护林面积（x₅）、水保固沙林面积（x₆）、用材林面积（x₇）、薪炭林面积（x₈）、经济林面积（x₉）、河谷低地草场面积（x₁₀）、沙地草场面积（x₁₁）、低山丘陵草场面积（x₁₂）、人工草场面积（x₁₃）、可利用水面积（x₁₄）、小畜头数（x₁₅）、大畜头数（x₁₆）及生猪头数（x₁₇）等，根据土地利用限制条件及所期望的几个目标值，可建立22个未标等级的“达成向量”。这些“达成向量”共分四类：

（1）第一类资源条件

①各类土地利用面积之和小于可利用面积：

式中：——第一个约束条件中的正偏差变量；

　　　——第一个约束条件中的负偏差变量。

②耕地面积小于多宜类土地面积：

③林地面积小于多宜类土地加上双宜类土地面积：

④放牧草场面积小于多宜类土地加上双宜类土地面积，再加上单宜类（宜牧）土地面积：

⑤本旗宜牧土地均为沙草地：

⑥各种土地利用所需要的劳力小于全旗劳力总数：

⑦用材林和经济林的宜种面积有一定限制：

⑧养鱼水面面积不超过可利用水面面积：

（2）生态要求

①森林面积不小于达到生态良性循环的最低数：

②防护林面积不小于生态要求的最低数：

③水保固沙林面积不少于生态规划数：

④牲畜总数不大于各种用地载畜量总和：

（3）社会需求

①粮豆总产不小于规划总产量：②油料总产不小于规划总产量：

③木材年采伐量不小于规划数：

④年采伐薪炭量不小于全旗农业人口所需能源的50%：

⑤大、小牲畜保持一定比例：

⑥生猪头数有一定保证：

⑦为解决过冬牲畜饲料，饲料地面积和牲畜头数保持一定比例：

⑧蔬菜、瓜果等其他作物面积保持一定数：

（4）经济要求

①年总投资不超过规划投资数：

②年总收＋人不小于规划数量：

根据以上“达成向量”及各组向量对偏差的要求，建立的模型目标函数为：

目标规划本身不能算出各个“达成向量”的决策变量系数，这些系数需要通过遥感、区划等途径获得数据，并预测2000年有关数据，采用数据平滑后的灰色系统模型群预测及精度估计，把灰色预测与数理统计有机结合起来，这就不仅能预测未来数据，还能估计预测精度。

设x⁰为原始数据列

将原始数据累加生成新的数据序列，考虑用下述形式的微分方程模型：

记系数向量a＝［a，u］^T，用最小二乘法对a求解：

从而建立起GM（1，1）的时间函数预测模型：

X^（1）（t＋1）＝（X^（0）（1）－U/a）e^－at＋U/a

当上述模型不能符合后险差的准确度要求时，还可以通过残差预测对模型进行补充修正。

建模采用的原始数据列X^（0）（t）的数目t≥4，假定收集到X^（0）有n个数，那么相连4个及4个以上的数据组一共可建立n²－n·（n＋3）/2个模型，从计算预测精度可比性要求出发，包含原始数据列中最后一位数的组合n－3，则可建立由n－3个子模型组成的灰色预测模型群。

以建立的第一个模型：

同理以建立的第n－3个模型。

以这些模型分别预测未来时刻t（t＝1，…，T）的值第t时刻预测值：

第t时刻预测精度为：

其中，

式中：a——信度，通常取0.05；

　　　n－4——自由度；

　　　——t分布可靠性指标。

本模型的22个“达成向量”分为5个目标等级。

资源条件是发展生产的前提，只有在绝对满足这一条件的情况下，才能考虑其他，列为第一级。一级目标向量共有8个，根据向量要求，目标函数为：

生态条件为第二级目标。主要包括森林覆被率、林地种类基本要求及草场载畜量等4个向量。

社会需要是发展生产的目的，列为第三级目标。包括农、林、牧产品等三个方面要求，共8个向量。

从经济效益来看，把投资作为第四级向量：

总产值作为第五级向量，则

以上述模型经计算机运算，获得奈曼旗土地利用最佳结构为：

耕地：151.1万亩（现有耕地178.20万亩），其中粮食作物102.27万亩，经济作物35.33万亩，饲料地10万亩，其他作物3万亩。

林地：470万亩（现有林地221万亩），其中防护林45万亩，水保防沙林240万亩，用材林129万亩，薪炭林30万亩，经济林26万亩。

草场：528.9万亩（现有草场603.7万亩）。

此外，养鱼水面10.8万亩，小畜50万头，大畜12.5万头，生猪30万头，在此基础上还对奈曼旗至2000年的粮食亩产及人均粮食、食油、木材、牲畜等占有量也作出了预测。

四、几点结论

1.为了使“三北”防护林遥感调查与专题制图建立在更高精度的基础上，首先必须有较好的信息源。针对“三北”防护林遥感调查的要求，通过对TM各波段信息特征的分析，并经过分段拉伸的TM_4（R）、TM_3（G）、TM_5（B）三波段合成为最优。通过实际处理，为“三北”防护林遥感提供了信息丰富、各种地物差异明显的图像，获得很好的解译效果。

2.对于“三北”防护林的重点地区，需要更大比例尺解译制图，可以通过计算机进行几何校正、空间卷积等处理，经C－4500光机扫描输出，可满足“三北”地区1∶5万专题制图的要求，甚至放大至1∶2.5万仍有较好的目视效果。与此同时，我们还提出了获得1∶5万～1∶2.5万高质量TM图像的关键技术方案。

3.NOAA卫星地面覆盖周期短，AVHRR资料通过比值运算，提取绿度信息和K－L变换处理再经密度分割及绿度指数与其他波段合成伪彩色处理，为“三北”防护林区的生态环境及防护林效益宏观分析和动态分析提供重要手段。

4.在遥感解译的基础上，进行土地资源评价机助自动评判及机助制图新方法的试验研究，为今后“三北”地区遥感动态监测及生态效益分析和资源评价提供了配套方法，为今后建立信息系统提供有利条件。

5.为了更好地研究“三北”防护林地区生产发展、生态平衡、合理地利用土地，我们把现状研究与今后规划发展预测研究结合起来，应用遥感数据对奈曼旗进行多目标规划研究，提出土地利用最佳结构，为当地生产发展和科学规划服务，并把遥感研究引向深入。

注：这篇文章是国家“七五”攻关“三北防护林区遥感调查”项目中“内蒙古草原防护林区”课题的技术总结，该课题获国家教委1993年科技进步一等奖。

【注释】

[1]参加本文资料收集工作的有：北京大学毛赞猷、崔海亭；内蒙古大学雍世鹏、李博、杨劼；华东师范大学陶然、宋浩昆、潘庆宁等，谨此致谢。