多光谱数据预处理

（一）多光谱图像的去噪声处理

实际上，噪声是一种高频信号，是以其灰度级与周围像元明显不同造成。不同图像产生噪声的强度及频率不同，因此，不可能有一种适合于所有图像的噪声消除算法。为此，我们选择了新疆阿尔金山地区一景有噪声的 TM数据，并设计了专门的程序来检测消除斑点噪声。具体算法是将一个像元的灰度值与周围像元进行比较，差值明显高于周围邻接像元，该像元就被定义为噪声，它的灰度值用邻近像元灰度值来替代或取周围与它距离最小的数个像元的算术平均值，且这数个像元灰度值不存在灰度突变。用这种方法可以消除单像元的孤立噪声点，且不会改变图像的信息量和清晰度（彩图1-1b）。对于面积分布较大的斑点噪声处理，可选择现有的几种主要方法进行，但利弊不同。其中，采用低通滤波可以有地效消除噪声，但同时也减少了其他有用的高频信息，为此对低通图像再进行高频滤波可以恢复部分高频信息。使用常规图像处理软件中的自适应滤波方法，如 Frost自适应滤波（彩图1-1c）、Lee自适应滤波、Gamma Map滤波等也有一定的效果，但图像的清晰度和信息量都有不同程度的降低。

（二）多光谱图像去条带处理

在TM图像上，条带的主要表现形式为：①图像数据行有若干相邻像元及整行丢失，图像上呈现白色条带。我们消除条带的具体做法是用“整行替代”算法，即用相邻行直接替代或用相邻两行的算术平均值去替代条带行（彩图1-2b）；②图像数据行有若干相邻像元及整行与周围行明显不同，其消除条带做法除了用“整行替代”的方法去消除之外，还可以采利用傅里叶变换提取水平条带，如在阿尔金山干旱裸露地区的实验中也取得了一定的效果（彩图1-2c）。具体方法是建立分解图像的傅里叶谱，穿过傅氏谱的垂直线表示水平条带，去掉这种垂直分量，然后再进行逆变换就可消除条带。

范正国

（地质矿产部航空物探遥感中心，北京　100083）

航空物探测量结果图件中的“条带”（下文中不加引号）现象，广泛地存在于各参量成果图件和转换处理图件之中，尤其是航空伽马能谱测量结果图件之中。图件中的条带，不但影响了图件的美观，更重要的是歪曲了地质构造的特征，给解释工作带来不便。因此，去条带的方法技术研究，一直受到国内外同行的广泛关注。已公开的方法主要有波纹滤波法（王懋基等，1991）、图像复原技术（张玉君等，1990）、大气氡日变观测及改正（水恩海等，1987）、相关分析法（A.A.Green,1987）、大气本底估值改正（Grasty,1986）、上测探头法（Foote,1968）、大面积水体上实测本底改正法（Darnley、Grasty,1970）、大气灰尘辐射改正法（Barson,1973）等。

条带的成因是相当复杂的，概略地说可包括两个方面，即不一致的测量条件（气象、温度、湿度、时间、大气放射线背景、飞行高度等）和不恰当的数据处理。尽管条带的成因是多方面的，但有一个共同的特征，即沿测线方向属长波异常，垂直测线方向属矩形波异常（阶梯状）。由于任何地球物理异常都具有连续渐变的特征，不存在阶梯状跃变的特点，故可以设计出特定的滤波器来去掉条带干扰。

一、方法原理

条带，实质是航空物探图件中出现的，沿测线方向延伸并具有一定宽度的条带状的场值升高或降低区段。其特点是，沿测线方向具有长波长特征，垂直测线方向具有阶梯状特征。换言之，条带属于区域性异常。由于区域地球物理场具有连续渐变的特点，不存在阶梯状跃变的特征。因此，从理论上讲，可以设计出特定的滤波器来去掉条带，同时又最大限度地保留区域地球物理场的固有特征。

（一）方向滤波法

由于条带是沿测线方向分布的长波长异常，因此，从理论上讲，可以通过方向滤波把条带异常分离出来。比较异常在空间域与频率域表现的主要特征，可以得出如下认识：空间域对水平叠加异常具有较大的分辨率，而频率域对垂直叠加异常具有较大的分辨率。因此，去条带的方向滤波法宜在空间域中进行。

设实测数据P（x、y）由真实数据f（x,y）和条带数据t（y）两部分组成，即

航空物探遥感论文集

式中，x为测线方向，y为垂直测线方向。用一维滤波算子H（x）作用于P（x,y），则有

航空物探遥感论文集

再用一维滤波算子H（y）作用于H（x）P（x,y），则有

航空物探遥感论文集

式（2）减式（3），并令

S（x,y）=H（x）f（x,y）-H（y）H（x）f（x,y）-H（y）t（y）

则有

航空物探遥感论文集

通过人机联作或统计分析的办法，可以找到一个合适的H（y）使S（x,y）→0，于是有

航空物探遥感论文集

把式（5）代入式（1），并移项，则得到数据改正的计算公式为

航空物探遥感论文集

（二）水平差商法

如前所述，条带是沿测线方向展布的形同条形码。因此，条带异常可以通过增强、计算机识别及人机联作等方法予以提取和剔除。

1.条带异常类型

从平面上看，条带与条形码类似。但条带数据与有用数据混合在一起，目前仍未有合适的滤波器能将二者彻底地分离开。为此，通过研究条带在垂直测线的横剖面（简称横剖面，下同）上的表现形式，进而提出识别、分离条带的方法。

在横剖面上，条带表现为阶梯状，沿横剖面求水平差商，则可得单脉冲异常、双脉冲异常、正负脉冲异常、多脉冲异常等，如图1所示。

图1　条带异常理想模式示意

1—无条带的原始场曲线；2—有条带的原始场曲线；3—无条带的水平差商曲线；4—有条带的水平差商曲线

2.条带异常的增强

实际工作中，除条带的水平差商能引起上述脉冲异常外，各种误差也能引起这类脉冲异常。为此，有必要增强源于条带的脉冲异常，压制误差引起的脉冲异常。

首先，需要沿测线方向对原始数据进行滤波，以减少随机误差的干扰。其次，横剖面应设置在平静场区，以减少有用异常的干扰。再其次，对水平差商做指数变换或幂变换，以增强条带异常（因为条带异常水平差商脉冲的幅值通常大于偶然误差产生的脉冲异常的幅值）。图2表示了一水平差商经幂变换前后的剖面曲线。变换前剖面曲线（1）存在明显的跳动弱异常；变换后剖面曲线（2）仅有条带异常，跳跃弱异常基本不存在。

3.条带异常的识别

通过上述增强及压制处理后，条带异常通常表现为强度大的陡直尖峰异常，而其它成因异常则为弱小的尖峰异常或宽缓异常，如图2所示。从该图上可以看出，目视方法能很好地区分条带异常与非条带异常。同样地，也可以设计适当的位变滤波器来识别条带异常。

图2　条带异常增强前后对比示意

1—增强前；2—增强后

4.条带异常的剔除

设条带异常两侧的水平差商值为△P1和△P2，则可用这两点的值对条带异常范围内的点做线性内插。待所有条带异常都内插完成后，再沿横剖面做积分，便得到了剔除条带后的区域场。为便于理解，举例说明如下。

设有一横剖面数据序列为｛—60，—50，—40,50,60，—30,0,10,20｝，其一阶向前差商△P为｛*,10,10,90,10，—90,30,10,10｝，剔除条带异常后的差商序列△f为｛*,10,10,10,10,10,10,10,10｝（*为无值），则剔除条带后的横剖面数据序列f为｛—60，—50，—40，—30，—20，—10,0,10,20｝。以上过程可用图3表示之。

以上实例表明，该方法不是通过常规滤波方法把条带碾平。因此，该方法能较好地保持真实数据的分布特征。图4为一条改正前后的真实横剖面，该图进一步表明了该方法的有效性。

图3　水平差商法去条带过程示意

a—原始横剖面Pb—原始水平差商△Pc—改正后水平差商△fd—改正后横剖面

图4　康古尔塔格地区横剖面曲线对比示意

1—原始剖面曲线；2—改正后的剖面曲线；3—改正前的水平差商曲线；4—改正后的水平差商曲线

5.剔除原始数据中的条带数据

设有M条横剖面P1（l=1,2，…，M），其改正后的数据为F1，则可求出区内第i线第j点的改正值为

航空物探遥感论文集

于是，第i线第j点改正后的值为

航空物探遥感论文集

式中pij、fij分别为第i线第j点改正前后的数值。

二、应用实例

针对上述方法研制了基于PC计算机的应用软件，并用其处理了青海柴达木盆地、新疆康古尔塔格地区航空伽马能谱测量数据以及甘肃潮水－雅布赖盆地、新疆伊犁地区航磁垂向一阶导数数据，效果良好。下面介绍在康古尔塔格地区的应用情况。

（一）康古尔塔格地区航空物探及地质简况

1990年航空物探遥感中心在新疆康古尔塔格地区完成了比例尺为1:2.5万、面积为6300km2的航空物探（磁、电、伽马能谱）综合测量，测线方向为南北，获得了质量较高的原始数据。尽管如此，在航空伽马能谱平面等值线图（图5）上仍存在明显的呈南北向展布的条带现象。该区主要地质构造格架（图6）是：北面以康古尔塔格深大断裂带为主体，其两侧次级断裂基本走向为北西、北东和近东西向；南面以雅满苏大断裂带为主体，其两侧次级断裂呈北西、北东向。区内放射性元素含量较高的地质体主要为中－酸性侵入岩，其长轴方向与本区断裂构造展布方向一致。另外，从航磁成果图件上亦可看出，本区不存在南北向的重大地质构造体。据此，航空伽马能谱图中的条带为非地质成因，应消除。

（二）去条带处理

限于篇幅，本文仅介绍消除航空伽马能谱钾含量数据中条带干扰的过程。

图5　康古尔塔格地区航空伽马能谱钾含量原数据等值线

图6　康古尔塔格地区地质构造略图

1—中酸性侵入岩；2—秋格明塔什—黄山韧性剪切带；3—阿奇山—雅满苏岛弧带；4—吐哈坳陷；5—断裂

1.用方向滤波法去条带过程

根据该区局部地质体南北方向长度（或宽度）可达到10km和条带异常一般均贯通南北，经反复试验，表明当选择沿测线方向滤波的滤波器半径为15km时，效果较好。

沿垂直测线方向滤波时，试验了不同滤波半径R情况下的去条带效果。当R为1km时，可去掉大部分条带（一般宽度较小的条带），同时，基本上不使区域场发生畸变。当逐渐增大R值时，去掉的条带越来越多，与此同时，区域场的畸变现象也逐步显现出来。当R达到9km时，图面上不再有条带，但区域场也出现明显的畸变。

2.用水平差商法去条带过程

首先，沿测线方向做滤波半径为0.7km的滤波，压制了噪声等可能对识别条带异常引起的干扰；然后将滤波后的数据做成图像，显示在计算机屏幕上；在没有明显异常处垂直条带设置横线切取剖面数据，将剖面数据求水平微商后用专用软件显示在计算机屏幕上（图4），再通过目视方法标出条带位置，计算机便自动做出校正；最后用图像显示、检查和修正结果。重复上述过程几十次，得到了满意的结果。

图7　康古尔塔格地区航空伽马能谱钾含量去条带后数据等值线

3.综合应用水平差商法与方向滤波法去条带过程

首先，重复 *** 作水平差商法去条带过程十余次，直到图像中不复存在宽度较大、强度较大的条带异常。然后，对上述结果沿测线方向做滤波半径为15km的滤波，求得区域异常。再沿垂直测线方向做滤波半径为0.8km的滤波，则可较好地消除条带异常。最后，用上述结果校正原始数据，得到校正后的数据，如图7所示。

对比图5、图6和图7可以发现，校正后数据较好地保存了原数据的基本特征，条带干扰不复存在，从而为解释人员提供质量更为可靠的原始图件。

三、结束语

为了较好地把航空物探数据图示出来，在成图之前，应消除包括条带在内的各种干扰。文中给出的去条带方法，经四个工区使用表明，①在较好地消除条带干扰的同时，较好地保存了原始基本特征。②水平差商法和方向滤波法各有利弊，前者保真度高但繁琐，后者便捷但保真度偏低。两者结合使用，用前者去掉宽度及幅度都较大的条带，用后者去掉宽度及幅度都较小的条带，其效果较好。

与国内外现有方法相比，有如下优点：①该方法不仅可以用于处理航空伽马能谱测量数据，也可以用于处理航磁数据；②由于该方法是对背景场进行修正，故不会造成异常畸变，这对保持原始数据的基本特征是很有必要的；③所研制软件是采用人机对话方式集成，并用彩色图像显示中间结果，从而使得去条带过程简单、明了，而且专家可以随时控制去条带的进程。

参考文献

1.熊光楚著.磁（重力）异常的变换及滤波技术.北京：冶金工业出版社，1990

2.王懋基等著.航空物探解释方法及应用.北京：地质出版社，1992

3.Green A A著.王红等译，利用道间关系修正航空伽马辐射数据.航测与遥感，1989,（4）:92～98

A NEW TECHNIQUE FOR ELIMINATING THE STRIPE INTERFERENCE IN AEROGEOPHYSICAL DATA

Fan Zhengguo

（Aerogeophysical Survey and Remote-Sensing Center,Beijing 100083）

Abstract

The stripe interference in aerogeophysical maps is a difficult problem which has not yet been satisfactorily solved.The present paper puts forward two new methods for eliminating stripe interference in aerogeophysical data.Based on the gradual change character of the regional geophysical field,the long wavelength feature of the stripe interference along the measuring line,and the step mutation and short wavelength feature of this interference vertical to the measuring line, this technique uses the intensification technique for stripe margins and computer（visual）recognition technique to extract stripe interference data and employs a special filter to eliminate such interference. Practice shows that this technique can not only eliminate the stripe interference but also quite satisfactorily retain the basic characteristics of the original data.

构建相山地区地学空间数据库是在对各类原始数据或图件资料进行整理、编辑、处理的基础上，将各类数据或图形进行按空间位置整合的过程。其工作流程见图 2.1。

图2.1 相山地区多源地学空间数据库构建流程

2.2.1 资料收集

相山地区有 40 多年的铀矿勘查和研究历史，积累了大量地质生产或科学研究资料。笔者收集的面上的资料包括原始的离散数据如航空放射性伽玛能谱数据、航磁数据、山地重力测量数据、ETM 数据，而地面高精度磁测资料仅收集到文字报告和图件。上述各类数据均可达到制作 1∶50000 图件的要求。地质图采用 1995 年核工业 270 研究所等单位共同实施完成的 “相山火山岩型富大铀矿找矿模式及攻深方法技术研究”项目的 1∶50000附图采用的 1∶50000 地形图的情况见表 2.1。

2.2.2 图层划分

GIS 数据库既要存储和管理属性数据和空间数据，又要存储和管理空间拓扑关系数据。数据层原理: 大多数 GIS 都是将数据按照逻辑类型分成不同的数据层进行组织，即按空间数据逻辑或专业属性分为各种逻辑数据类型或专业数据层。相山地区数字化地质图包括地理要素和地质要素两大部分，共设置 9 个图层，每一图层 (包括点、线或多边形) 自动创建与之相对应的属性表。

表2.1 采用的地形图情况一览表

注: 坐标系均为 1954 年北京坐标系，1956 年黄海高程系，等高距为 10 m。

(1) 水系图层 (L6XS01) : 包括双线河流、单线河流、水库或水塘。

(2) 交通及居民地图层 (L6XS02) : 包括公路和主要自然村及名称。

(3) 地形等高线图层 (L6XS03) : 包括地形等高线及高程和山峰高程点。

(4) 盖层图层 (D6XS04) : 包括第四系 (Q) 和上白垩统南雄组 (K2n) 及其厚度和主要岩性。

(5) 火山岩系图层 (L6XS05) : 包括下白垩统打鼓顶组 (K1d) 、鹅湖岭组 (K1e) 及各种浅成- 超浅成侵入体 (次火山岩体) 的分布和主要岩性特征。

(6) 基底图层 (L6XS06) : 含下三叠统安源组 (T3a) 、震旦系 (Z) 、燕山早期花岗岩 (γ5) 、加里东期花岗岩 (γ3) 。

(7) 构造图层 (L6XS07) : 相山地区褶皱构造不发育，构造图层主要包括实测的和遥感影像解译的线性断裂或环形构造。

(8) 矿产图层 (L6XS08) : 包括大、中、小型铀矿床和矿点。

(9) 图框及图幅基本信息图层 (L6XS09) : 数字化地质图的总体描述，内容包括图框、角点坐标、涉及的 1∶500000 标准图幅编号、调查单位及出版年代等。

图层名编码结构如下:

相山铀矿田多源地学信息示范应用

2.2.3 图形输入

图形输入或称图形数字化，是将图形信息数据化，转变成按一定数据结构及类型组成的数字化图形。MapGIS 提供智能扫描矢量化和数字化两种输入方式。本次采用扫描矢量化输入，按点、线参数表事先设定缺省参数，分别将地形底图和地质底图扫描成栅格图像的 TIF 文件，按照图层划分原则，在计算机内分层进行矢量化。线型、花纹、色标、符号等均按 《数字化地质图图层及属性文件格式》行业标准执行。

对于已建立的图层，按点、线、多边形分别编辑修改，结合地质图、地形图及相关地质报告，采集添加有关属性数据，用以表示各图层点、线、多边形的特征。拓扑处理前先将多边形的地质界线校正到标准图框内进行修改，去掉与当前图层区域边界无关的线或点。对于图幅边部不封闭的区域，采用图框线作为多边形的边界线，使图幅内的多边形均成为封闭的多边形。拓扑处理后进行图形数据与属性数据挂接。

在 MapGIS 实用服务子系统误差校正模块中，将数字化地图校正到统一的大地坐标系统中。图形数据库采用高斯-克吕格 (6 度带) 投影系统，椭球参数: 北京54/克拉索夫斯基。

MapGIS 数据文件交换功能使系统内部的矢量图层很容易实现 Shape 和 Coverage 等文件格式的转换。在图形处理模块将上述各图层转成 Shape 文件格式。

2.2.4 离散数据网格化

在收集的原始资料中，除 1∶50000 地形图和地质图之外，航空放射性伽玛能谱数据(包括原始的和去条带处理后的数据) 、航磁数据、山地重力测量数据都是离散的二维表格数据。用 GeoExpl 网格化。GeoExpl 数据处理与分析系统提供了多种网格化计算的数学方法，本次选用克立格插值方法，网格间距 15 m。重力和航磁数据网格化后，进行不同方向或不同深度的延拓处理。所有网格化数据均采用了与上述图形数据相同的地图投影和坐标系统。

2.2.5 网格化数据影像化

MapGIS 网格化文件格式为 grd，可直接被 Erdas Imagine 读取，GeoExpl 网格化文件包括重磁处理反演后的网格化文件可转换成 Surfer.grd 后，被 Erdas Imagine 读取。然后将上述网格化数据一一转成 img 影像数据格式。

2.2.6 DEM 生成

地形等高线 (L6XS03) 文件在 MapGIS 空间分析子系统 DEM 分析模块中，生成 DEM栅格化文件: L6XS03.grd，再转成 img 格式，文件名改为: XSDEM。

经过上述程序形成的各类矢量或栅格数据，在 ArcView 平台建立 “相山数据库”工程文件，将上述各 Shape 图形和 img 影像文件一一添加到该工程文件中。该工程文件即为相山地区矢量、栅格一体化地学空间数据库。该数据库，一可以对这类地学空间信息实现由 GIS 支持的图层管理，二可以视需要不断进行数字—图形—图像的转换，三可以将多源地学信息进行叠合和融合，以实现多源地学信息的深化应用和分析，为实现相山地区铀资源数字勘查奠定基础。

---