遥感图像压缩

2024-06-07

遥感图像压缩（共8篇）

遥感图像压缩篇1

摘要：压缩感知理论,以远低于奈奎斯特采样定理要求的速率对图像进行采样,可利用图像的部分信息重构原始图像,有效地减轻图像处理的计算复杂度,降低对硬件的要求。文中提出了一种基于压缩感知的遥感图像融合方法。在压缩域对多光谱和全光谱遥感图像进行了融合实验,并与传统的融合方法进行了比较,实验结果表明,文中方法在遥感图像融合上有着良好的性能。

关键词：压缩感知,图像融合,小波分解,图像重构

随着传感器系统的快速发展,如何从传感器数据提取科学信息成为重点。在许多情况下,一个传感器所提供的信息无法反映现实情况,因此,在过去的几年里多传感器融合引起了业内关注。图像融合是多传感器融合的一个分支,是指由两个或两个以上的图像结合成一个拥有超过输入图像的任何一幅的信息的融合图像的过程。随着图像的数量增加以及传感器的发展,对于过程中图像的获取和硬件的要求越来越高。近几年来提出的压缩感知理论(Compressed Sensing,CS)[1,2]标志着图像处理领域一个崭新时代的来临。在CS理论中,图像的采样和压缩同时进行,可以用远低于奈奎斯特采样定理要求的速率对图像采样,而不会丢失图像的任何有用信息,达到完全重构原始图像的目的。该方法可以减少信号采样率、存储量、功耗、图像处理及相关研究领域的计算复杂性,有效地降低了对传感器硬件的要求。

基于CS图像融合的自然方法就是在图像融合后采用CS方法处理。但为减少计算复杂度和节省存储空间,更好的方法是将融合规则应用于经过压缩后产生的测量值。Tao Wang等[3]提出了对测量值小波变换后进行一种简单取大融合的规则,Xiaoyan Luo等[4]提出了对测量值基于熵的加权平均的融合规则,Juanjuan Han等[5]对文献[3]中的方法进行了改进,并与文献[4]中的方法进行比较,得到了更为优越的结果。

文中在提出了一种基于压缩感知理论的图像融合框架,在对测量值的融合中,引入了梯度的方法,融合规则采用小波分解低频取全色系数,高频取梯度大的系数,有效提高了图像融合的效果,并与已有算法进行了比较。实验结果表明,文中算法有着较好的优越性。

1 压缩感知理论

1.1 压缩传感

压缩传感的核心是线性测量过程,设x(n)为传统采样得到的稀疏信号,长度为N,通过压缩传感得到y(m),长度为M,M<N,它们的关系为y=Φx,其中Φ称为测量矩阵,大小为M×N。如果x(n)不是稀疏信号,可在正交稀疏变换下通过系数s(k)表示,记为y=Θs,其中Θ=ΦΨ为M×N矩阵。

显然,无法直接从y(m)恢复出x(n),因为上述方程组中未知数个数N超过方程个数M,但可通过求解下列优化问题得到x(n)在Ψ变换域内的稀疏形式可压缩传感的关键是通过构造合理的测量矩阵,采用非适应的重建算法对信号进行恢复。从以上可以看出,M值越小,数据压缩也就越大,但M越小,数据就越难恢复,所以他们之间有个平衡量,如果要得到99%以上的正确恢复率,三者之间应满足M≥4K(ln N+4.6)。同样,对于测量矩阵,如果Φ满足RIP特性,即(1-δK)x22≤Φx22≤(1+δK)x22,0<δK<1,就可以恢复出原始信号。

1.2 重构算法

文中采用OMP重构算法[6,7,8],其本质思想是:以贪婪迭代的方式选择Φ的列,使得在每次迭代中所选择的列与当前的冗余向量最大程度地相关,从测量向量中减去相关部分并反复迭代,直到迭代次数达到稀疏度m,强制迭代停止。具体算法步骤如下:

(1)初始化冗余向量r0=y,索引集合Λ0=Ø,迭代计数t=1。

(3)令Λt=Λt-1∪{λt}。

(4)计算{φλ∶Λt}张成空间的正交投影Pt。

(5)计算新的近似at和冗余rt∶at=Pty,rt=y-at。

(6)t=t+1,如果t<m,返回步骤(2)。

(7)获得的估计sλ^在索引Λm位置的元非零,且在该位置的值就是对应的分解系数:。

2 基于压缩感知理论的遥感图像融合

2.1 基于CS的遥感图像融合方法

文中基于CS理论的图像融合框架主要思想[9]是,首先分别对输入的源图像进行压缩采样,并在感知域上进行图像融合,然后将融合观测传输给接收端,接收端通过OMP算法进行重构,得到输入源图像的融合图像。该框架融合过程如图1所示。

上述融合方法的主要包括信号调制、融合和图像重构3大步骤:

步骤1信号调制。将输入的两幅图像x1和x2分别利用同一欠采样观测矩阵进行观测,得到观测图像y1和y2,也就是将图像低维投影到CS感知域上。

步骤2融合。在CS感知域上,将步骤1中获得的观测图像y1和y2依据一定的融合规则进行融合,得到融合观测图像y。

步骤3图像重构。利用步骤2得到的融合观测图像y,基于CS理论进行图像重构,得到最终融合图像x^。

2.2 基于CS的遥感图像融合算法

在图像处理中,多尺度小波分解对于信号表示有着显著优势[10,11]。小波分解的基本思想是将任意函数f表示为小波的迭加,这种函数f的小波迭加表示将f分解为不同尺度级,在每个尺度级,f又在这一尺度级对应的分辨下被分解,尺度级对应着频率,频率越高,对应的分辨率越高,换而言之,将函数f分解成不同频率通道成分,并将每一频率通道成分又按相位进行了分解——频率越高者,相位划分越细,反之则越疏。

文中采用二维小波分解算法[12],其定义:令Vj2(j∈Z)为L2(R2)的可分离多分辨逼近空间,φ(x,y)=φ(x)(y)为对应的二维尺度函数,ψ(x)为与尺度函数φ(x)对应的一维标准正交小波。

设待分析图像f=f(x,y)∈Vj2,则二维逼近图像Ajf(x,y)可分解为

用Hr和Hc分别表示对阵列{Ck,l},(k,l)∈Z2的行和列作用的算子H,用Gr和Gc分别表示算子G作用于该阵列的行和列上,二维Mallat重构算法为

文中采用二维小波分解得到逼近系数CA,细节水平系数CH,细节垂直系数CV,细节对角系数CD。文中的融合算法为:逼近系数取全色图像,细节系数分别在不同方向求其梯度,最后分别取梯度最大值为其细节系数。通过二维小波逆变换得到融合后的观测向量,再通过OPM重构算法得到融合结果,并将融合结果逆稀疏化后得到融合图像。

3 基于压缩感知的遥感图像融合实验

实验由信息熵、平均梯度、光谱扭曲度及光谱相关系数来反映实验结果。融合图像信息熵的大小反映了融合图像空间细节信息量的多少。定义为

平均梯度反映了图像中微小细节反差表达的能力和纹理特征的变化。定义为

光谱扭曲度是融合图像灰度与源多光谱图像灰度差值图像的平均值。定义为

光谱相关系数表示融合图像与源图像光谱特征的相关程度。定义为

4 实验结果

源图像为256×256的全色图像和多光谱图像,采用小波稀疏化,观测矩阵为哈达玛矩阵,重构采用OMP算法,由于文献[5]中已经与文献[3～4]方法进行了比较,结果表明文献[5]方法性能更好。文献[5]采用一维小波分解的方法,文中方法采取二维小波分解,实验中与文献[5]中的1d方法进行了比较。实验结果如图4所示。

从图4中可以看出,文中融合方法有着良好的融合结果,而且从客观数据中也能看出其优越性。

从表1可以看出,在信息熵、平均梯度、光谱扭曲度、光谱相关系数等方面,文中方法都优于文献[5]的融合结果。可见文中提出的方法在提高空间信息和保留多光谱信息之间找到了最佳平衡点,即以最小的空间信息量下降为代价,换取了更多的多光谱信息保留。

5 结束语

压缩感知理论对图像处理领域具有重要意义。文中提出了一种基于CS理论的遥感图像融合方法,并通过实验参数说明了该方法的优越性。但是,算法仍存在一定的不足,文中并没有对CS域的重构算法及测量矩阵进行进一步的研究,这将是后期进一步探索的方向。

遥感图像压缩篇2

学院理学院班级地信131 学号姓名

编写日期：1

2015.5

▶▶作业a

1.LS8_C_20140613_022505_000000_118039_GEOTIFF_L4

2.L5118_39_19860531 ProductDescription用记事本打开，读取头文件，并填写相关信息与相应位置即可

3.L5118-39-19960103

4.L7118039_20050815 直接打开以_mtl为后缀的文件，该文件中包含了遥感影像的所有波段 5.LM*** 直接打开波段，然后波段合成即可

6.s5kj297_289_10m

7.WORLDVIEW-052606622010_01

▶▶作业b

在ENVI中将landsat的4景影像和SPOT-5的1景的影像打开，并联动连接查看同一区域

link displays是根据象元位置来连接的，geographic link是通过地理坐标位置来连接的。

由上图可知，将遥感影像联动时亦可实现不同影像同一区域的快速检索，但是我们也可以看到，由于受到各方面因素的影像并不能特别精确的指在同一地方。

▶▶作业c

1.WORLDVIEW-2影像保存为jpg和TIF格式的4-3-2波段合成的假彩色图像。可用同样的方法将SPOT-5影像保存为jpg和tif格式的4-3-2波段合成的假彩色图像 2.为landsat的5景影像附上波段的波长，并根据波长用landsat 5的7-4-3波段，保存为jpg和tif格式影像

为波长复制后，导入影像文件各波段显示差异前后对比

转换为JPG格式后可以用看图软件直接打开

▶▶作业d

需要对影像进行裁剪，裁剪的基本步骤如下：

1.L5118_39_19860531裁剪前后对比

2.L5118-39-19960103裁剪前后对比

3.L7118039_20050815裁剪前后对比

4.LS8_C_20140613_022505_000000_118039_GEOTIFF_L4裁剪前后对比

▶▶

作业e

将剪裁影像，重采样成10m，重采样的操作主要如下

1.L5118_39_19860531重采样前后对比

2.L5118-39-19960103重采样前后对比

3.L7118039_20050815重采样前后对比

遥感图像压缩篇3

1 遥感图像的特点

遥感图像分辨率较高,目标物体所占的像素较少,图像中所含的高频信息和纹理信息较多,空间频率的局部变化较快,图像的不连续性较强。

为了从空间相关性来分析遥感图像数据之间的特征,文中采用图像的自相关函数r(l,k)来衡量

r(l,k)=∫∫f(x,y)f(x+l,x+k)dxdy (1)

其中,f(x,y)为图像的灰度函数。对上式进行归一化和离散化处理,有

$R (l, k) = \frac{\sum_{x = 1}^{Μ} \sum_{y = 1}^{Ν} [f (x, y) - u_{f}] \times [f (x + l, y + k) - u_{f}]}{\sum_{x = 1}^{Μ} \sum_{y = 1}^{Ν} [f (x, y) - u_{f}]^{2}} (2)$

其中,图像的大小为M×N,uf为图像亮度的平均值。即

$u_{f} = \frac{1}{Μ \times Ν} \sum_{x = 1}^{Μ} \sum_{y = 1}^{Ν} f (x, y) (3)$

当l=k=1时,即表示相邻两个像素之间的自相关函数值R(1,1),称为自相关系数。

在对遥感图像进行3层小波分解后,除LL子带集中了主要的低频信息,其它所有高频子带近似遵循零均值Gamma 分布。分析表明,除极少数纹理较强的局部区域之外,绝大多数纹理不太强的区域近似服从广义Gaussian 分布,特别是平坦区域,均值大多约为0。遥感图像的这些特性为采用JPEG2000中的算法提高图像压缩效率提供了理论依据。

2 遥感图像高频子带分级量化方案

遥感图像进行JPEG2000整数小波变换后,编码系统需要对分量中的每个拼接块(Tile-component)中的系数进行量化。量化的关键是根据变换后图像的系数特征、重构图像质量要求等因素设计合理的量化步长。其中高频子带影响较小,只作粗量化编码;低频子带对重建过程影响较大,可进行精细量化编码。

2.1 基于纹理特征的高频子带分级量化方案设计

本方案的具体步骤如下[1,2]:

(1) 对图像进行n层小波分解,将遥感图像分解成3n+1个子带,将子带进行层次划分,依次分成LLn,LHn-HLn,HHn-LHn-1,…,HH2-LH1-HL1和HH1,如图1所示;

(2) 将LLn子带直接送入编码器进行EBCOT 编码和算术编码;

(3) 依次将除LLn子带之外的高频子带小波系数划分为4×4的块,设定阈值T1>T0>0。将方差<T0的子块划分为平坦区,方差>T0且<T1的子块划分为次纹理区,方差>T1的子块划分为纹理区;

(4) 将LHn-HLn层中的小波系数分别进行分级量化。取量化步长1≤Δ $_{1}^{(n)}$ <Δ $_{2}^{(n)}$ Δ $_{3}^{(n)}$ 分别对子带中各区域进行量化。对块内方差<T $_{0}^{(n)}$ 的平坦区,采用量化步长Δ $_{3}^{(n)}$ 对小波系数进行粗量化;对块内方差>T $_{1}^{(n)}$ 且<T $_{1}^{(n)}$ 的次纹理区域,采用量化步长Δ $_{2}^{(n)}$ 对小波系数进行量化。对块内方差>T $_{1}^{(n)}$ 的纹理区域,采用量化步长Δ $_{1}^{(n)}$ 对块内小波系数进行精细量化;

(5) 将HHn-LHn-1-HLn-1层中的小波系数进行分级量化。依次对HHn,LHn-1,HLn-1进行区域量化处理,其量化步长分别为Δ $_{1}^{(n - 1)}$ <Δ $_{2}^{(n - 1)}$ <Δ $_{3}^{(n - 1)}$ 且Δ $_{1}^{(n)}$ ≤Δ $_{1}^{(n - 1)}$ 。设定其阈值为T $_{1}^{(n - 1)}$ >T $_{0}^{(n - 1)}$ >0,具体的处理过程与步骤(3)相同;

(6) 将HH1子带中的小波系数进行量化。由于HH1子带所含的能量极少,可以直接取量化步长Δ(0)(Δ $_{1}^{(1)}$ ≤Δ(0))对其小波系数进行量化处理,或对该子带做清零处理;

(7) 在完成上述步骤后,将获得的系数进行舍入处理整数化,获得新的小波系数。再将其划分为编码块(Code-block)送入Ebcot和算术编码器中进行编码;

(8) 将各高频子带对应的量化步长依次写入相应的拼接块(tile)标头的QCD或QCC字段中,各子带的阈值T0和T1也依次写入相应拼接块的标头com字段中,最终形成JPEG2000码流格式。

在解码器端,由于采用本方法压缩后,码流组织为JPEG2000码流结构,可以将其直接进行JPEG2000解码,获取恢复图像[3]。

2.2 关于阈值和量化步长的选取

(1) 阈值T0,T1的选取。

阈值的选取可以有两条途径:一是通过对遥感图像进行实验可以获得各个阈值与压缩系统性能指标峰值信噪比和压缩比的关系,由用户根据需要交互地给出,或由系统给出一个较优的缺省值;另一种方法是通过对小波分解子图进行统计分析后自适应获得[2,4]。

(2) 量化步长的选取。

量化步长的选取可以根据子带的统计特性和最小失真,计算出量化步长的最优值;实际应用中可先选取初始量化器,根据实际图像压缩过程,验证各量化步长与压缩性能指标PSNR和CR(峰值信噪比和压缩比)的关系,经过多次调整量化步长最终确定较优值。

3 仿真结果及分析

选取卫星遥感图像中截取的512×512,8 bit灰度图像作为测试图像。其中图2(a)为航拍遥感图像,图2(b)为卫星遥感图像,如图2所示。

小波变换滤波器为Daubechies(5/3)滤波器组,进行3级离散小波变换。在分级量化方法中还附加条件如下

T1=2×T0 (4)

${\begin{cases} Δ_{1}^{(3)} = 1 \\ Δ_{2}^{(i + 1)} = Δ^{(i)} \\ Δ_{j + 1}^{(i)} = Δ_{j}^{(i)} \times 2 \end{cases}$

其中,i∈{0,1,…,3},j∈{1,2,3} (5)

表1给出了上述两幅遥感图像采用JPEG2000无损压缩后的测试数据。由于无损压缩,图像可以完全重构,PSNR趋近于无穷大。所以采用压缩重构时间、压缩比为图像压缩的评价标准,具体数据,如表1所示。

表2给出了以上两幅遥感图像采用高频子带分级量化后的测试数据。采用分级量化方案对遥感图像进行压缩,由于高频子带采用>1的量化步长,势必产生一定的失真;另外,分级量化后的参数整数化处理也引入了部分失真。因此,在评价此方案的压缩性能时采用压缩和重构时间、峰值信噪比和压缩比为评价标准。

由表1和表2可知,与JPEG2000无损压缩方法相比,分级量化方法在图像压缩性能上有了明显提高,不同图像压缩比提高了20%左右;尽管引入了部分高频子带系数的失真,峰值信噪比有所下降,但是采用文中方案重建的图像质量均可达到40 dB,人眼很难分辨出图像的失真。因此,本方案能满足遥感图像压缩的要求。

采用分级量化方法对遥感图像进行高保真压缩和重构,运算时间均有所增加。由于在编码过程中,将高频子带分块,计算每个子块的方差并对其进行量化;解码过程中,需要判断高频子带子块的纹理属性,也通过计算方差并与阈值进行比较。因此,此方法在一定程度上增加了算法的复杂度。

4 结束语

文中在JPEG2000无损压缩的框架上,增加了一种遥感图像高频子带分级量化方案,实验表明其实现简单,压缩效果良好。可见,该方法在遥感图像高保真压缩方面具有良好的效果。

摘要：介绍了JPEG2000静态图像压缩标准的基本情况;对遥感图像的相关性以及小波变换域统计特性进行了分析,得出遥感图像的特点。针对其特点,设计了基于JPEG2000无损压缩的高频子带分级量化高保真压缩方案。实验结果表明,此方案在对遥感图像进行压缩时取得了较好的效果。

关键词：遥感图像,图像压缩,小波变换

参考文献

[1]Chao H,Fisher P.An Approach of Integer Reversible Wave-let Transform for Image Compression[J].The Wavelet Di-gest,1996,5(8):19-38.

[2]万登峰.基于JPEG2000遥感图像压缩的应用研究[D].上海:上海交通大学,2006.

[3]王达生.ISO/IEC/JTC1/SC29/WG11第53次会议将召开[J].计算机辅助设计与图形学学报,2000,12(5):374-380.

遥感图像压缩篇4

关键词：遥感图像,小波包,父冲突,扩展方向树,SPIHT

1 引言

遥感,是20世纪60年代兴起并迅速发展起来的一种信息高技术,已在城市规划、环境保护、军事侦察等领域得到了广泛的应用[1]。遥感图像的特征是数据量庞大,图像尺寸高达30000×30000像素,谱段包含数百个,因此,要保证遥感数据能够在有限带宽的信道中传输,就必须进行压缩。

与常规图像相比,遥感图像的纹理更多、分辨率更高,空间相关性很小,采用普通的压缩方法可能造成纹理、边缘的损失,因此,如何在提高压缩比的同时又能较好地恢复细节也就成为了遥感压缩的重点。基于此,我们利用小波包优良的高频分析能力[2],提出了一种结合SPIHT的小波包编码算法—扩展方向树小波包压缩算法(Extended SOT Wavelet Packet Coding)。该方法借鉴SPIHT算法中的零树结构[3,4],通过重新定义其方向树,即扩展方向树(Extended SOT),解决零树结构中的“父冲突”(parental conflict)[7,8,9]问题,最终通过SPIHT算法实现整个编解码。

2 小波包变换

小波包变换[2]是小波变换的推广,它通过对小波变换的高频部分作进一步分解,使信号频带划分更加精细,并且能够自适应地选择小波包基,以达到对原始信号的最优逼近。小波包定义如下所示:

称由式(1)所定义的函数族{un(n=0,1,2···,n∈Z+)}为由φ生成的小波包。其中u0(t)=φ(t),u1(t)=ψ(t),φ(t)和ψ(t)分别为尺度函数和小波函数。

同时,将un进行二进伸缩和平移可以得到:

定义由式(2)组成的函数库{uj,n,k(t),j,k∈Z,n∈Z+}为由φ(t)导出的小波包库。由式(2)不难看出,小波包库由许多小波包组成,而不同小波包又具有不同性质,反映不同的信号特性,小波包必然具有分解形式多样化的特点。因此,本文采用Shannon熵[5,9]作为代价函数来得到小波包变换的最优基。

3 扩展方向树

传统小波编码如SPIHT算法,成功地利用了小波变换的零树结构和子带间相关性,通过方向树有效地表示了有效值映射,其方向树定义如图1所示。

图中每个子带用黑色正方形表示,其边长与尺度成反比。除最低频子带LLM和最高频的三个子带外,每个系数在它的同方向相邻更高频子带的相同位置上有四个子女,最高频子带没有子女,同方向相邻子带间尺度满足二进关系。

然而,由于小波包有可能对高频子带作进一步分解,相邻子带间差异可能不再满足二进关系,其父子关系会更复杂,子系数可能出现在比父系数尺度更粗糙的子带中。这样就出现了如下两类“父冲突”:

1)一个子系数对应多个父系数。2)不同节点的父系数对应同一个子带的子系数。

以三层小波包分解为例(图2(a)所示),R表示最低频子带LL3,T1,T2,T3表示水平,垂直,对角三个不同方向树的根节点,每个方向在下一精细尺度的子带都是其子节点。这样,可以定义如下四种父子关系:

a)子节点位于父节点同方向相同位置的更粗糙尺度子带。如图2(b),树T2:14与15、16、17、18所示。

b)子节点位于父节点同方向相同位置的相同尺度子带。如图2(b),树T2:14与19、20、21所示。

c)子节点位于父节点同方向相同位置的下一个精细尺度子带。如图2(b),树T3:7与11所示。

d)子节点位于父节点同方向相同位置的更精细尺度子带。如图2(b),树T1:1、2、3、4与5所示。

可以看出,A、D两种父子关系与小波变换父子关系不同,产生了“父冲突”问题,我们称A为第一类“父冲突”,D为第二类“父冲突”。对于前者,文献[8]通过将产生冲突的粗糙尺度的子节点上移,使之与更粗糙或者相同尺度下的节点产生父子关系而解决冲突。而对于第二类“父冲突”,由于四个父节点都对应于同一个子节点,我们将其拆分,产生冲突的子节点只与其中一个父节点产生父子关系。综合考虑C、D两种父子关系,一个父节点就不是仅仅对应一个子节点,一个父系数也不是仅仅对应4个子系数,而是应该对应4N,N表示父子之间的尺度差异。因此本文采用的方向树不同于SPIHT算法,我们称为扩展方向树,在实际编码时需对子节点个数和子孙数目做相应调整。图3所示解决“父冲突”后的方向树图。

4 零树假设验证

对于一个给定的小波包分解,其方向树的定义如上所述。而零树结构是建立在一个假设的基础上的:在小波分解系数中,对于给定的阈值,如果父系数是不重要的,那么很有可能子孙系数也是不重要的[3,4]。正是利用这种假设,才能够有效表示“无效值映射”,也就等价于有效地表示了有效值映射,这就是零树结构的基础。因此,本文编码算法是否成功,首先就取决于零树假设在小波包分解中是否成立。我们引入小波分解零树结构检验的经验统计方法,该方法同样适用于小波包分解。

我们采取两种方式统计变换系数的幅值:1)按照子带所属频带,以频率递增的方向画出所有系数幅值;2)对不同的方向树,从根节点开始,按照先父节点后子节点的顺序分别画出水平,垂直以及对角方向的系数幅值。以一幅2048×2048的遥感图像为例,小波及小波包分解后子带扫描顺序分别如图4(a)、图4(b)所示,扫描按照子带标号递增方向进行,子带内部则以行为单位,按照从上到下的顺序扫描。图5则表示在当前小波包分解下所生成的扩展方向树。

图6和图7分别表示小波分解和小波包分解后按频率增加顺序和方向树顺序画出的幅度值。在图6(b)和图7(b)中,小波分解和小波包分解都能够按照方向树将系数分为三个部分,而且系数基本呈现由大到小逐渐变化的趋势。对于小波包分解,在每一个树内,还存在若干子树(图8),每一个子树内部同样呈现父节点系数幅度较大,子节点系数幅度较小的趋势。换言之,对同一个树,父系数总是比子系数大,如果父系数不重要,子系数也很有可能不重要,因此,小波包扩展方向树满足零树假设。

5 实验

为验证ESOT-WPC算法的准确性,我们选择3幅2048×2048×8bit的遥感图像进行实验,并与SPIHT算法比较,三幅图像分别命名为RS1、RS2和RS3。小波包分解选用双正交Db9/7小波,5层小波分解,并通过代价函数的计算确定最优基,解决方向树中的“父冲突”问题,最优基结构如图9所示。SPIHT算法同样选用双正交Db9/7小波,5层小波分解。

图10、图11和图12分别是RS1、RS2和RS3在0.5bpp压缩率下,使用SPIHT算法和ESOT-WPC算法压缩后得到的图像。不难看出,图10(a)与图10(b)和图10(c)在视觉观察上基本无异,这说明如果图像冗余度较高,包含大量平坦区域,本文算法和S PIHT算法都能得到较好的效果;对于RS2和RS3,图11(a)中边缘在图11(b)中虽有所恢复,但已经发生扭曲,而在图11(c)中则基本保持了边缘的线性。图12(a)白色柱状物之间的斜向纹理在图12(c)中基本得到恢复,但在图12(b)中则是模糊一片。表1是三幅图像使用SPIHT算法和ESOT-WPC算法压缩后的PSNR(峰值信噪比)的对比,对于细节变化很少的RS1图像,由于图像中包含大量的低频信息,过多的高频分解反而会降低图像的压缩效率,因此本文算法的P SNR低于SPIHT。而对于富含纹理和细节的图像RS2和RS3,本文算法在无论在PSNR还是在视觉效果上都明显优于SPIHT算法,更好的恢复了遥感图像的细节信息。

6 结论

本文针对普通压缩算法容易造成遥感图像高频信息丢失的问题,提出一种结合SPIHT的小波包编码算法—ESOT-WPC。借鉴SPIHT算法中的零树结构,通过重新定义子带之间的父子关系,即扩展方向树,解决其“父冲突”问题,并成功实现编解码。实验表明,该算法在PSNR和视觉效果上均优于SPIHT算法,是一种有效的遥感图像压缩方法。但目前算法计算复杂度较高,需进一步改善算法性能,降低复杂度。

参考文献

[1]彭望琭.遥感概论[M].北京:高等教育出版社,2002.PENG Wang-lu.The Overview of remote sensing[M].Beijing:Higher Education Press,2002.

[2]Ramchandran Kannan,Vetterli Martin.Best wavelet packet bases in a rate distortion sense[J].IEEE Trans.Image Processing,1993,2(4):160-175.

[3]Said Amir,Pearlman William A..A new fast and efficient image codec based on set portioning in hierarchical trees[J].IEEE Trans.on CASVT,1996,6(3):243-250.

[4]Shapiro Jerome M.Embedded image coding using zero-trees of wavelet coefficients[J].IEEE Transactions on Signal Process-ing,1993,41(12):3445-3462.

[5]Coifman Ronald R,Wickerhauser Mladen Victor.Entropy based algorithm for best basis selection[J].IEEE Trans.Information Theory,1992,38(2):713-718.

[6]Xiong Zixiang,Ramchandran Kannan,Orchard Michael I.Wavelet packet image coding using space-frequency quantization[J].IEEE Trans.Image Processing,1998,7(6):892-898.

[7]Xiong Zixiang,Ramchandran Kannan,Orchard Michael I d.Space-frequency quantization for wavelet image coding[J].IEEE Trans.Image Processing,1997,6(5):677-693.

[8]Rajpoot Nasir M,Wilson Roland G,Meyer Francois G,et al.Adaptive wavelet packet basis selection for zerotree image coding[J].IEEE Trans.Image Processing,2003,12(12):1460-1472.

遥感图像压缩篇5

在航天应用中,线阵CCD应用广泛,因为它具有较面阵CCD更多的优点,它的结构简单,单排感光单元数可以做得很多,在同等测量精度下,它的测量范围可以做得较大,而且能够在低照度条件下工作。多通道输出的线阵CCD能够在相对较短的时间内实现图像数据的转移,相对地提高了CCD的感光时间,但是输出的图像数据首先要进行数据整合才能得到一幅完整的图像。

图像压缩在各个领域中应用广泛,通过压缩可以减少图像的存储空间,提高传输线路的利用率,而无损图像压缩在航天、航空侦查以及医学摄像领域中应用广泛,较之于有损图像压缩方法,它能够将原始图像无失真地还原出来。JPEG-LS无损图像压缩算法是一种算法难度低、压缩效率较高的无损压缩算法,它的算法难度低便于硬件实现。采用FPGA来实现该算法,可以实现流水线的工作模式,提高算法的执行效率。

1 系统的总体结构

整个系统的工作原理如图1所示,CCD传感器经过视频处理器处理后的数据分为4个通道输出(如图中的A、B、C、D),在双口静态存储器的协助下,将四个通道输出端的图像数据整合成一副完整的图像数据,紧接着输入到JPEG-LS无损压缩模块[1]中进行图像压缩,源图像以预定模式的输入给编码器,在经过上下文建模的判别后决定采用哪种方式进行压缩:一种是常规(Regular)模式,即对图像进行预测,得到预测残差,然后对预测残差进行哥伦布(Golomb)编码[2];另一种是游程(Run)模式,即对图像直接进行游程编码, 得到的压缩图像通过一个输出端口输出。

2 CCD数据的整合

2.1 四通道的线阵CCD工作原理

四通道输出的12 000像元的线阵CCD,如图2所示,A、C通道为奇像元通道,B、D通道为偶像元通道。A通道输出1～5 999奇数像元的数据,B通道输出2～6 000偶像元通道数据,C通道输出6 001～11 999奇数像元数据,D通道输出6 002～ 12 000偶数像元数据。C、D通道输出的数据是排在后面的像元数据先输出,排在前面的像元数据后输出的,这与A、B通道数据输出的顺序相反。

2.2 CCD数据整合的FPGA实现

本文的FPGA选用Xilinx公司的Spartan3的XC3S200芯片,双口静态存储器选用IDT公司的IDT70V27,它的大小为32 K×16 bit。

图 3所示是系统仿真的时序图,CCDCLK是像元转移时钟(为5 MHz),FRAME是帧转移始终,DataA-DataD是四个通道输出的数据,CLKRam是驱动双口RAM的时钟,WAddr是写RAM的地址,RAddr是从RAM读数据的地址。因为双口RAM能同时进行读写操作,为了保证数据能够不间断地处理,在双口RAM中开辟两块区域,一块用于将CCD数据写入RAM中, 另一块是从RAM中读取数据。在向RAM中写入数据时,要求在一个像元转移时钟内,写入四个通道的数据,并且写A,B的数据从第一个地址递增,写C,D通道的数据是从第12 000个地址开始并且依次递减。从RAM中读取数据时, 初始地址为相对于写操作的另一块地址区域的起始地址并且依次递增,这样就读出了一幅完整的图像,用于图像压缩。

3 FPGA实现JPEG-LS图像压缩

3.1 JPEG-LS压缩方法

JPEG-LS压缩方法有两种工作模式:一种是无损压缩模式;另一种是有损压缩模式。本文用的是它的无损压缩模式,在无损压缩工作模式下,JPEG-LS压缩算法的简化框图如图4所示,源图像按预定的扫描模式依次输入给编码器,编码器采用如图中左边的模板的来预测x的值,并通过a,b,c,d的梯度变化来决定编码器采用常规模式编码还是采用游程编码模式,采用常规编码模式时,预测器通过两个部分来预测x的值,一部分是固定预测,即仅根据a,b,c,d的值来预测出一个值,然后通过自适应修正器来对这个预测值进行修正,用x处的真实图像值减去预测值,得到预测误差,通过上下文建模器,对预测误差进行调整,使得预测误差在自己规定的范围内,然后通过哥伦布编码器来进行编码输出压缩后的图像。在游程模式中,编码器直接跳过预测和误差编码,从x处开始计数其值和a的值相同的一连串像素,即游程的长度。

3.2 JPEG-LS无损压缩实现

3.2.1基于FPGA实现JPEG-LS的结构与原理

实现JPEG-LS算法的流程为:

1)初始化。初始化365个预设的上下文计数值A,B,C和N,初始化游程编码的表以及索引,降抽样值x指向图像的第一个像素;

2)计算模板的梯度值。g1=d-b;g2=b-c;g3=c-b; 判断g1=g2=g3=0 ,若等于零则转向游程编码(第八步),若不等于零则为常规模式,继续下一步;

3)量化g1,g2,g3的梯度值为q1,q2,q3,并将其映射到范围为[1]其中的一种上下文条件下Q;

4)根据以下公式计算固定预测值 $\hat{x}$ MED,然后自适应修正器对固定预测值进行修正得到 $\hat{x}$ :

$\hat{x}_{Μ E D} ≜ {\begin{cases} \min (a, b) i f c \geq \max (a, b) \\ \max (a, b) i f c \leq \min (a, b) \\ a + b - c o t h e r w i s e \end{cases}} (1)$

5)计算预测误差,即将真实像素值减去预测值 $\hat{x}$ 得到 $\bar{ε}$ ,并将 $\bar{ε}$ 映射到预定范围内;

6)计算哥伦布编码的参数k,将预测误差进行哥伦布编码;

7)更新A,B,C,N的值,并返回第一步进行下次循环;

8)游程扫描,对游程长度进行编码。

在FPGA内采用流水线[3]的方法来实现算法。

将算法分成五个模块,各个模块的工作复杂度几乎相同,流水线工作模式如图5所示。在JPEG-LS算法描述里,变量表的更新是在最后一部完成的,但是下次操作要用到这些变量表,所以在FPGA中,将旧的变量A,B,C,N的值保存,便于后面的操作,而同时更新A,B,C,D的值,为了下一次计算,从而保证了流水线工作模式的完成。在每个模块中尽量采用并行算法的模式,如“计算预测误差,保存旧的A,B…0,更新A,B…的值”这个模块,可以将计算预测误差,保存和更新A,B…两个步骤并行运算,减少每一级的运算时间。

3.2.2 设计实现

设计采用大量的片内RAM资源,首先是对A,B,C,N变量表的存储,使用片内Distributed Ram资源产生4个730字节的变量表,将梯度量化值q1,q2,q3,归一化为Q值,将其转换成索引存储器的地址,实现对A[Q],B[Q],C[Q],N[Q]的索引。当预测误差(Errval)计算出来之后使用更新A,B…模块更新这个变量表里的值。其次,在计算梯度g1,g2,g3时,模板需要两行图像数据,所以要对前一行数据进行存储,在初始化时由于还没有图像数据,将该行的数据初始化为零,在计算梯度的同时,存储读取的图像数据为正在处理的行,若当前行处理完毕,则将该行的地址作为下一行处理的前一行数据的地址。这里设计采用片内的Block RAM 资源设计一个24 000 bytes的双口RAM。这样设计需要的RAM资源占片内资源的88%。

量化g1,g2,g3为q1,q2,q3这个过程如下(这里用C语言描述,Verilog语言也类似),其中T3,T2,T1,NEAR为梯度阈值:

//JPGE-LS标准 //FPGA中实现qi相应的值

If (gi <=-T3) qi=-4; 0;

Else if(gi <=-T2) qi=-3; 1;

Else if(gi <=-T1) qi=-2; 2;

Else if(gi <-NEAR) qi=-1; 3;

Else if(gi <=NEAR) qi=0; 4;

Else if(gi < T1) qi=1; 5;

Else if(gi < T2) qi=2; 6;

Else if(gi < T3) qi=3; 7;

Else qi=4; 8;

在JPEG-LS标准中量化出来的值按照0对称,在FPGA中为了能够更简单的实现查找A,B,C,N中的表,将量化的qi做了修改,在归一化到Q时采用以下公式:

$Q = q_{1} + q_{2} \times 9 + q_{3} \times 9^{2} (2)$

A,B,C,N的变量表按照所需要的Q的升序来进行安排时,Q就可以作为索引这些变量的存储器地址线了,在算法中当q1为负(即SIGN=-1)时,q1,q2,q3都取反,算法要求有365种上下文(Contexts)情况。这里需要730个地址来索引,其中一半重复的,这样增加了存储数量来减少复杂度和提高效率。

使用公式(1)计算出x的固定预测值(Fixed prediction)px,通过以下公式对px进行调整。其中MAXVAL是算法限定的最大预测值。C[Q]即在Q处的C的值:

$p_{x} = {\begin{cases} p_{x} + S Ι Ν G \times C [Q] \\ Μ A X V A L (p_{x} > Μ A X V A L) \\ 0 (p_{x} < 0) \end{cases}} (3)$

得到预测值之后便可以使用真实的像素值减去预测值,得到预测误差Errval,得到预测误差之后,将A[Q],B[Q],C[Q],N[Q]的值备份,存为A1[Q],B1[Q],C1[Q],N1[Q]。然后更新计算A[Q],B[Q],C[Q],N[Q]的值,并将这四个存储区在Q地址处的值更新。

计算Golomb编码的k值,k值的计算是将N1[Q]值右移k位,当N1[Q]的值小于A[Q]时,这个k值即Golomb编码所需要的k。将预测误差Errval重新映射,使得误差为非负数。

进行Golomb编码,以2k作为Golomb编码的参数。表1所示为k=2时,Golomb编码值。

在FPGA中使用Verilog语言实现q=(Errval-1)〉〉k;r=x-(q〈〈 k)-1;编码值为((((1〈〈q)-1)〈〈1)〈〈k)|r。

游程编码在FPGA中的实现的Verilog HDL语言如下,以RUNval+RUNcnt的编码就是游程编码:

4 仿真与结果

整个设计使用Verilog HDL语言[4]实现,在Xilinx的ISE开发平台下完成,综合后占用8 400个逻辑单元,220 k的RAM资源。结合Modelsim对多幅图像进行了仿真,仿真结果表面,对一个像元输入直到预测误差输出(在系统时钟为20 MHz的条件下)延时约为260 ns,压缩一幅12 000像元的图像为1.2 ms。

5 结论

根据FPGA能并行运算的优势,将JPEG-LS无损压缩的结构,按照流水线的模式进行分割,为设计硬件实现图像的无损压缩提供了较好的解决方案,无损压缩能在使图像不失真的条件下,减少了图像的传输量,为图像的传输系统减轻了负担。

参考文献

[1]Marcelo Weinberger,Gadiel Seroussi,Guillermo Sapiro.TheLOCO-I lossless image compression algorithm:principles andstandardization into JPEG-LS[J].IEEE Trans Image Proc,2000,9(8):1310-1322.

[2]Golomb S W.Run length encoding[J].IEEE Transactions onInformation Theory,1966,12:399-401.

[3]Rafael C Gonzalez,Richard E Woods.Digital Image Proces-sion[M].Second Edition.BEIJING:Publishing House ofElectronics Industry,2002:326-400.

遥感图像分类方法比较研究篇6

1.遥感图像分类原理简介

遥感图像分类的理论依据是:遥感图像中的同类地物在相同的条件下 (纹理、地形、光照及植被覆盖等) , 应具有相同或相似的光谱信息特征和空间信息特征, 从而表现出同类地物的某种内在相似性, 即同类地物像素的特征向量将集群在同一特征空间区域; 而不同的地物其光谱信息和空间信息特征不同, 将集群在不同的特征空间区域。由于地物的成分、性质、分布情况的复杂性和成像条件, 以及一个像素或瞬时视场里往往有两种或多种地物的情况, 即混合像素, 使得同类地物的特征向量不尽相同, 而且使得不同地物类型的特征向量之间的差别不都是截然相反的。

遥感图像分类方法可分为非监督分类和监督分类。非监督分类是在没有先验类别知识的情况下, 根据图像本身的统计特征及自然点群的分布情况划分地物类别的分类处理, 它并不需要具体地物的已知知识, 可更好地获取目标数据内在的分布规律, 主要方法有:K均值分类 (K-means) 、ISODATA方法。

监督分类是一种常用的精度较高的统计判决分类, 它是在已知类别的训练场地上提取各类训练样本, 通过选择特征变量、确定判别函数或判别规则, 把图像中的各个像元划归到各个给定类的分类方法。其主要步骤包括:选择特征波段;选择训练区;选择或构造训练分类器;对分类精度进行评价。主要方法有:平行六面体分类 (Parallelepiped) 、最小距离分类 (Minimum Distance) 、马氏距离分类 (Mahalanobis Distance) 、最大似然分类 (Maximum Likelihood) 、人工神经网络分类法等。

2.案例分析

本文以呼和浩特市地区7个波段组成的TM影像为原始数据, 以波段5、4、3模拟真彩色RGB合成, 选取凉城县岱海区域做分类实验, 以确定该区域土地覆盖应采取的分类方法。

2.1遥感影像的目视解译。

在利用ENVI进行监督分类之前, 首先在ArcGIS中进行遥感影像的矢量化操作, 利用目视解译画出各类地物的边界, 并统计其面积, 以此评定各种分类方法的分类精度。

2.2遥感图像的平行六面体法、最小距离法、马氏距离法、最大似然法分类。

采用ENVI软件监督分类模块进行平行六面体法、最小距离法、马氏距离法、最大似然法分类, 分类时训练区相同, 以比较不同分类方法的精度。分类后对原始分类图像进行聚类、过滤、去除后处理工作, 并采用总体分类精度、Kappa系数 (结果略) 评定分类精度, 发现最大似然法分类精度最高。图1~图4为四种方法原始与后处理分类图像。

2.3遥感图像的人工神经网络分类。

2.3.1在ENVI主菜单中选择Classification→Supervised→Neural Net Classification, 弹出Classification Input File对话框, 选择分类影像shilianghua.tif。点击“OK”打开Neural Net Parame- ters对话框, 如图5所示。

2.3.2参数设置 :Select Classes from Regions:单击Select All Items按钮, 选择全部的训练样本。

Activation:选择活化函数。备选函数为对数 (Logistic) 和双曲线 (Hyperbolic) 函数, 此处选择对数函数;Training Threshold Contribution:输入训练贡献阈值 (0~1) 。该参数决定了与活化节点级别相关的内部权重的贡献量, 它用于调节节点内部权重的变化, 训练算法交互式地调整节点间的权重和节点阈值, 从而使输出层和响应误差达到最小。将该参数设置为0, 不会调整节点的内部权重, 适当调整节点的内部权重可以生成一幅较好的分类图像, 但是如果设置的权重太大, 对分类结果就会产生不良影响, 此值需多次试验选择;Training Rate:设置权重调节速度 (0~1) 。参数值越大则训练速度越快, 但也增加摆动或者使训练结果不收敛;Training Momentum: 输入一个0~1的值。该值大于0时, 在“Training Rate”文本框中键入较大值不会引起摆动, 该值越大, 训练的步幅越大, 该参数的作用是促使权重沿当前方向改变;Training RMS Exit Criteria: 指定RMS误差为何值时训练应该停止。RMS误差值在训练过程中将显示在图表中, 如图6所示。当该值小于输入值时, 即使没有达到迭代次数, 训练也会停止, 然后开始分类;Number of Hidden Layers:所用隐藏层的数量。要进行线性分类, 键入值为0 (没有隐藏层) , 不同的输入区域必须与一个单独的超平面线性分离。要进行非线性分类, 输入值应该大于或等于1, 当输入的区域并非线性分离或需要两个超平面才能区分类别时, 必须拥有至少一个隐藏层才能解决这个问题, 两个隐藏层用于区分输入空间, 空间中的不同要素不临近也不相连;Number of Training Iterations:输入用于训练的迭代次数。这里选默认值1000;Min Output Activation Threshold:输入一个最小输出活化阈值。如果被分类像元的活化值小于该阈值, 在输出的分类中, 该像元将被归入未分类 (unclassified) 中。

2.3.3输入文件输出路径 , 单击“OK”执行分类 , 分类结果如图7所示。

2.4五种遥感图像分类方法比较。

将目视解译、平行六面体、最小距离、马氏距离、最大似然及神经网络分类法统计的面积进行对比, 结果如图8所示。

将五种分类结果与目视解译的结果求差, 可看出神经网络分类结果与目视解译结果最为接近。差值图如图9所示。神经网络分类精度比最大似然法高, 且高于其他三种分类方法。

3.结语

如何提高遥感影像的分类精度, 一直是关注的热点。本文通过对五种监督分类方法的对比分析, 发现人工神经网络土地利用分类与实际情况最接近, 分类精度最高。最大似然法次之, 两者分类精度皆高于平行六面体、最小距离法、马氏距离法分类方法。

人工神经网络与传统分类方法相比, 有其明显优点。

一方面, 神经网络方法没有对数据分布特征的任何假设前提。它不考虑数据是正态分布或是不连续分布, 可以在特征空间上形成任意边界的决策面, 并在每次迭代过程中动态调节决策区域, 显示统计方法强大的稳定性和优越性。神经网络的监督分类和标准的统计监督分类非常相似。主要区别在于训练和分类在实现的细节上。神经网络的训练过程和最大似然法计算均值和方差类似。所不同的是, 最大似然法的均值和方差只计算一次, 而神经网络法采用迭代算法, 直到计算结果和实际结果的差异满足要求。

另一方面, 神经网络方法是非线性的, 同传统分类方法相比较, 它可以处理复杂的数据集, 识别精细的模式, 且能利用多源数据, 将潜在的信息提取出来。

遥感图像压缩篇7

关键词：遥感,几何校正,主成份分析,数据融合,专题图像,ENVI

0 引言

ENVI是一个完善的数字影像处理系统,它具有全面分析卫星和航空遥感影像的能力。它能在各种计算机操作平台上提供强大新颖的友好界面,显示和分析任何数据尺寸和类型的影像。由于采用了基本文件和基于波段的技术,ENVI能够处理整个影像文件或单一的某个波段。此外,ENVI自带了多种波谱处理工具,包括波谱提取、波谱库以及高光谱数据集的能力,这些高光谱数据集涵盖了AVIRIS、GERIS、GEOSCAN和Hy Map等数据。

1 几何校正

遥感成像时,由于飞行器姿态(侧滚、俯仰、偏航)、高度、速度、地球自转等因素而造成图像相对于地面目标而发生几何畸变,畸变表现为像元相对于地面目标实际位置发生挤压、扭曲、伸展和偏移等,针对几何畸变进行的误差校正就称为几何校正。

将数据参照水系图进行几何校正,该方法是:

(1)坐标转换从原始畸变图像空间中的像元位置出发,建立空间转换关系,确定每个像元在标准图像空间中的正确位置。

(2)重采样用标准图像空间中的像元点G位置反求其在原始畸变图像空间的共轭点F(X,Y),然后再利用某种方法确定这一共轭点的灰度值,并把共轭点的灰度值赋给标准图像空间对应点g。

重采样法能够保证校正空间中网格像元呈规则排列,因而是最常用的几何精校正方法之一。重采样有3种插值方法:最邻近法(nearest neighbor)、双线性法(bilinear)、三次卷积法(cubic convolution)。

最邻近法保留大量原始灰度值,没有经过平滑处理,对于区分植被类型、识别线性特征等有重要意义。双线性内插法比与最近邻法相比,计算量增加了,但提高了精度,改善了灰度不连续现象及线状特征的块状化现象。其缺点是这种方法对图像起到平滑作用,使图像变得模糊。由于这种方法计算量和精度适中,因而常常被采用。

在控制点的选取中选择道路的交叉点,河流弯曲或者分叉处,及固定地物边缘等。在这个过程中我选取了16个控制点,其中RMS Eerror分别为0.28(X)和0.46(Y),数学模型选用Polynomial方法对图像进行校正处理。

选取地面控制点(GCP,Ground Control Point)建立多项式,多项式纠正法的精度与地面控制点的精度、分布、数量及纠正范围有关;GCP的位置精度越高,则几何纠正的精度越高;GCP的个数不少于多项式的系数个数;适当增加GCP的个数,可以提高几何纠正的精度。GCP分布应尽可能在整幅图像内均匀分布,否则在GCP密集区精度较高,在GCP分布稀疏区出现较大误差;GCP应该选择明显地物点,如河流交汇点、桥梁、道路拐点等。

2 运用校正后的数据求NDVI并叠加

NDVI(Normalized Difference Vegetation Index)是一个普遍应用的植被指数,将多波谱数据变换成唯一的图像波段显示植被分布。NDVI值指示着像元中绿色植被的数量,较高的NDVI值预示着较多的绿色植被。

步骤:Transforms→NDVI(Vegetation Index),出现NDVI Calculation Input File窗口,选择输入文件,点击“OK”确定,通过点击“Input File Type”下拉菜单,用NDVI Calculation Parameters对话框,说明你已经输入的文件类型,用“Output Data Type”下拉菜单选择输出类型,选择输出到“File”或“Memory”,若选择输出到“File”,在标有“Enter Output Filename”的文本框里键入要输出的文件名;或用“Choose”按钮选择一个输出文件名,点击“OK”开始计算NDVI变换。先是得到NVDI图,然后在此基础上得到其叠加图。

3 数据融合

最常用的遥感影像融合方法有HIS(亮度(intensity)、色调(hue)、饱和度(saturation)变换、brovey方法、主成分变换(principal component transform,简称PC)、小波变换、基于亮度调节的平滑滤波(smoothing filter-based intendity modulation,简称SFIM)和Gram-Schmidt光谱锐化(Gram-Schmidt spectral sharpening,简称Gram),等。这些方法各有优缺点,根据不同需要选择不同的方法。

在PCI的Modelor实现8个波段(123457,61-62,8波段除外,因为8波段具有高分辩率,在HIS变换中需要再次用,故此处不用它来参与)PCA变换,HIS变换,变换后的图像进行图像增强,最后输出文件。

PCA变换提取8个波段的所有有用信息;用提取出的PC1,PC2,PC3当成R,G,B3个通道分别输出到RGB to HIS conversion模块中;在接下来的HIS to RGB conversion模块中用saturation和hue分别输出到saturation和hue中,用part8高分辩率图像作为intensity数据输入,输出后的数据就是经过HIS融合后的数据。

ETM+的8个波段波段主要的用途有:(1)0.45~0.52um蓝绿波段,用于水体穿透,土壤植被分辨;(2)0.52~0.60um绿色波段,用于植被分辨;(3)0.63~0.69um红色波段,处于叶绿素吸收区域,用于观测道路/裸露土壤/植被种类效果很好;(4)0.76~0.90um近红外波段,用于估算生物数量,尽管这个波段可以从植被中区分出水体,分辨潮湿土壤,但是对于道路辨认效果不如TM3;(5)1.55~1.75um中红外波段,这被认为是所有波段中最佳的一个,用于分辨道路、裸露土壤、水,它还能在不同植被之间有好的对比度,并且有较好的穿透大气、云雾的能力;(6)10.5~12.5um热红外波段,感应发出热辐射的目标。分辨率为60m;(7)2.08~2.35um中红外波段,对于岩石/矿物的分辨很有用,也可用于辨识植被覆盖和湿润土壤;(8)0.52~0.90um全色波段,得到的是黑白图象,分辨率为15m,用于增强分辨率,提供分辨能力。(ETM+增加了这个波段,TM没有)。

以上波段除6、8外分辨率均为30m。

选择Transform->Image Sharpening->HSV,以校正后的Band8图像为High Resolution Input File,融合校正后的Band 123457与校正后的Band8图像。

4 利用融合后的多光谱数据做主成分分析,生成PC1、2、3的彩色图像,并对结果进行解释

彩色增强就是根据人的视觉特点,将彩色用于图像增强之中,这是提高遥感图像目标识别精度的一种有效方法。彩色合成增强是将多波段黑白图像变换为彩色图像的增强处理技术。根据合成影像的彩色与实际景物自然彩色的关系,彩色合成分为真彩色合成和假彩色合成两种。真彩色图像上影像的颜色与地物颜色基本一致。利用数字技术合成真彩色图像时,是把红色波段的影像作为合成图像中的红色分量、把绿色波段的影像作为合成图像中的绿色分量、把蓝色波段的影像作为合成图像中的蓝色分量进行合成的结果。真彩合成。与肉眼所见接近;仅使用反射的可见光,受大气、云雾、阴影、散射的影响较大,通常对比度不高,感觉模糊(蓝色光散射严重);对于海岸区域研究特别有用,因为可见光可穿透水面,观察到海底。

以下是一幅专题图像:(包含专题图像、地名标注、图框、坐标、比例尺、图例等要素)

非监督分类中采用k-means分类将图像分为5类。

由植被分布专题图结合地面实况资料,广东地处我国岭南,具有高温多雨而有较明显的干、湿季之分的季风热带气候特点;这里属低平地上的地带性典型植被类型,植被及其生境的热带-亚热带过渡性质。

图像分类是对图像上每个像素按照亮度接近程度给出对应类别,以达到大致区分遥感图像中多种地物的目的。遥感图像的统计分类一般分为两种方法:监督分类和非监督分类。在此我用非监督分类的方法做此图。

5 总结

遥感影像作为地表地物光谱信息的一种载体出现,不仅扩大了人们的视觉领域,更是提高了人们认识世界、认识自然的能力,给生产科研带来深刻的革命性变化。利用ENVI4.2,对遥感数据及其它地学数据在不同层次上进行的多波段数据融合、不同图像数据融合、图像校正、图像分割、分类处理,其结果表明,利用ENVI的遥感图像处理技术在突出地质地貌特征信息方面有较好的效果。

通过本次实例分析,我从中学习到很多遥感图像处理的方法,但是还有以下不足:

实习非监督分类中训练区选取不太适当,所以分类的结果不太理想,在以后的分析过程中应该针对需要的类别选取特征性的训练区。分类中具体方法针对性不同,选择分类效果更好的方法。

参考文献

[1]梅安新,彭望路,秦其明.遥感导论[M].北京:高等教育出版社.

[2]宁书年.遥感图像的处理与应用[M].北京:地震出版社.

[3]赵英时.遥感应用分析原理与方法[M].北京:科学出版社.

[4]夏德深.现代图像处理技术与应用[M].南京:东南大学出版社.

遥感图像几何校正原理及效果分析篇8

在二十一世纪的今天, 遥感已揭开她神秘的外纱, 与多种应用学科相结合, 从理论遥感已发展出资源遥感、环境遥感、生态遥感、海洋遥感、地质遥感等领域, 遥感应用方面的教材和专著也层出不穷。本人多年从事遥感理论与应用方面的研究和教学工作, 发现遥感图像处理中的几何校正是一个易被忽视的研究内容, 常造成处理后遥感图像的质量下降。本文结合实例, 较深入地剖析几何校正原理, 较全面地指出提高几何校正精度的有效参数和途径。重视几何校正在遥感图像处理中的作用, 它不仅可以有效提高遥感图像的几何精度, 而且可以进行遥感图像与其它数据源如GIS图形的配准。

1 遥感图像几何校正原理分析

在分析遥感图像几何校正原理的同时, 有必要了解几何校正产生的背景以及与遥感图像几何精度的关系。

1.1 几何畸变与几何校正

几何畸变和几何校正是遥感理论的一对派生词, 几何校正是因几何畸变而产生, 是解决几何畸变的方法体系。

在遥感理论上, 将遥感平台位置和运动状态、地形起伏、地球表面曲率、大气折射等遥感系统内外因素影响造成的遥感图像几何位置上的变化统称为几何畸变, 也就是遥感图像在几何位置上与实际地面位置有差异。在图像上表现为像元行列分布不均匀, 像元大小与地面大小对应不准确等。针对不同因素产生的图像几何畸变, 发展出了多种多样的校正方法, 如图像比例尺的变化、中心位置的偏移以及旋转等畸变, 通过不同的线性变换进行校正。对于图像的不规则几何畸变, 通过多项式法、三角测量法、MQ模型、自动纠正算法等众多方法[1,2,3,4]来消除。

1.2 图像级别

遥感图像的几何精度级别与几何校正密切相关。一般遥感图像数据源分为3个级别, 第1个级别是指卫星下行数据经过格式化同步、按景分幅、格式重整等处理后得到的图像数据。该数据未进行任何几何校正, 几何畸变较大。第2个级别是经过系统几何校正的图像数据, 即利用卫星所提供的轨道和姿态等参数、以及地面系统中的有关处理参数对第1个级别遥感数据进行几何校正所获得的数据, 该工作通常是由遥感卫星地面接收站完成, 是研究人员获取的主要遥感数据源。第3个级别是经过几何精校正的图像数据, 指采用地面控制点对第2个级别中的系统几何校正模型进行修正, 使之更精确地描述卫星与地面位置之间的关系。该工作可以由遥感卫星地面接收站完成, 也可以由研究人员来实现。

对于第3级别的Landsat 7 ETM+遥感图像数据, 几何定位精度可达一个像元以内, 即几何误差为15米, 可以作为1:10万以下中小比例尺的矿产、资源、测量等研究方面的数据源。

1.3 几何校正原理分析

几何校正, 即图像级别中的几何精校正, 主要方法是多项式法。方法机理是通过若干控制点, 建立不同图像间的多项式空间变换和像元插值运算[5,6], 实现遥感图像与实际地理图件间的配准, 达到消减以及消除遥感图像的几何畸变, 是获取第3个级别遥感数据的重要方法。该方法适用于地面平坦, 不考虑高程信息, 或地面起伏较大而无高程信息, 以及传感器位置和姿态参数无法获取等情况。同时在遥感图像制图、不同类型或不同时相遥感图像间的几何配准和复合分析上也有重要应用。如人们习惯于正射投影地图, 可以采用以不同正射投影地图为无几何畸变图件或标准图, 对非正射投影图像进行几何校正, 得到正射投影的遥感图像专题图。

多项式几何校正机理实现的两大步:

(1) 图像坐标的空间变换

有几何畸变的遥感图像与没有几何畸变的遥感图像, 其对应像元的坐标是不一样的, 如图1右边为无几何畸变的图像像元分布图, 像元是均匀且等距分布;左边为有几何畸变的遥感图像像元分布图, 像元是非均匀且不等距的分布。为在有几何畸变的图像上获取无几何畸变的像元坐标, 需要进行两图像坐标系统的空间转换。

在数学方法上, 对于不同二维笛卡儿坐标系统间的空间转换, 通常采用的是二元n次多项式, 表达式如下:

${\begin{cases} x = \sum_{i = o}^{n} \sum_{j = 0}^{n - i} a_{i j} u_{i} v_{i} \\ y = \sum_{i = o}^{n} \sum_{j = 0}^{n - i} b_{i j} u_{i} v_{i} \end{cases}$

其中x, y为变换前图像坐标, u, v为变换后图像坐标, aij, bij为多项式系数, n=1, 2, 3, …。

二元n次多项式将不同坐标系统下的对应点坐标联系起来, (x, y) 和 (u, v) 分别对应不同坐标系统中的像元坐标。这是一种多项式数字模拟坐标变换的方法, 一旦有了该多项式, 就可以从一个坐标系统推算出另一个坐标系统中的对应点坐标。

如何获取和建立二元n次多项式, 即二元n次多项式系数中a和b的求解, 是几何校正成败的关键。数学上有一套完善的计算方法, 核心是通过已知若干存在于不同图像上的同名点坐标, 建立求解n次多项式系数的方程组, 采用最小二乘法, 得出二元n次多项式系数。

不同的二元n次多项式, 反映了几何畸变的遥感图像与无几何畸变的遥感图像间的像元坐标的对应关系, 其中哪种多项式是最佳的空间变换模拟式, 能达到图像间坐标的完全配准, 是需要考虑和分析的。

在二元n次多项式数字模拟中, 从提高几何校正精度的角度考虑, 需要兼顾的因素主要有引起几何畸变的原因和产生数学运算误差因素。归纳起来有三个方面的考虑因素:

一是多项式中n值的选择, n值与几何畸变的复杂程度密切相关。当n=1, 上述的坐标空间变换成为二元一次多项式, 可以进行线性的坐标变换, 解决比例尺、中心移动、歪斜等方面的几何畸变, 实用于第2级别以上的遥感数据。n值的不同选择, 可以得到不同的空间变换式, 当n≥2, 上述的坐标空间变换成为二元非线性多项式, 解决遥感器偏航、俯仰、滚动等因素引起的几何畸变。从理论上讲, n值越大, 越能校正复杂的几何畸变, 但计算量也相对要大。实际应用中n值通常取小于等于3。

二是控制点GCP (用于空间坐标变换的同名坐标点) 的选择, GCP的几何精度直接影响着多项式系数的求解误差大小。成熟的作法是:通过目视, 选择熟悉的、易分辩且精细较高的特征点 (如小水塘边缘、线状地物的交叉点、海岸线弯曲处等) , 且GCP分布在全图中要尽量均匀, 特征变化性大的地区选择多些, 图像边缘部分选些控制点, 使系数的求解尽可能准确。控制点精度的衡量尺度为RMS (Root Mean Square) 参数, 意为均方根, 以图像像素大小为单位, 表达式为:

$R Μ S = \sqrt{(x - x^{'})^{2} + (y - y^{'})^{2}}$

x, y为无几何畸变的图像控制点坐标, x′, y′为变换后图像控制点坐标。在ERMAPPER7.0或ENVI等遥感软件中, 对于一次线性变换, 当采集到4个控制点以上时, 软件系统就会自动推算控制点的变换值和RMS, 可以很好地辅助控制点的编辑。在实际应用中需要引起注意的是:随着控制点数目的增减, 多项式系数值也在变化, 每个控制点的RMS也在变化。当RMS值都小于等于1时, 控制点的精度控制在一个像素大小上, 几何校正效果较好。

最后是控制点GCP数目的确定, 从数学运算上来说, 一次多项式变换, 存在6个系数要计算, 需要GCP的最少数目是3。二次多项式变换, 有12个系数需要计算, GCP最少数目是6。n次多项式, GCP的最小数目为 (n+1) (n+2) /2。但在实际应用中, 采用最小GCP数目, 几何校正效果往往不好。所以在条件允许的条件下, GCP数目要远远大于最小数目, 可以是其的6倍。

(2) 图像像元灰度值重采样

经过上述图像像元坐标的空间变换, 得到对应于实际地面或无几何畸变的图像坐标, 图像上每个像元都有了无几何畸变的坐标值。随后需要做的是给每个像元赋亮度值。因为已知的图像数据是有几何畸变的像元亮度值, 并没有校正后的无几何畸变的像元亮度值。所以需要通过数学上的重采样方法如最近邻法、双向线性内插法和三次卷积内插法等计算出校正后像元位置的亮度值, 形成无几何畸变的遥感数据。在重采样方法中, 三次卷积内插法计算量虽最大, 但图像质量要好, 细节表现要清晰, 是许多遥感软件的首选方法。

2 广东省从化市Landsat 7 ETM+几何校正

2.1 数据准备

本文选择广东省从化市2003年11月9日Landsat 7 ETM+遥感图像数据和近期同地区1:50万水系GIS图件。Landsat 7 ETM+图像 (如图2所示) 含有8个波段, 1、2、3、4、5、7波段分别为可见光、近红外和短波红外波段, 波长依次为0.4787, 0.5610, 0.6614, 0.8346, 1.65, 2.208μm, 地面分辨率为30m。8波段为全色波段, 波长为0.5-0.9μm, 地面分辨率为15m。6波段为热红外波段, 波长为10.4—12.5μm, 地面分辨率为60m。

以1:50万水系GIS图件 (如图3所示) , 包括本地区近期1级至6级的水系分布, 为标准的无几何畸变图像, 对ETM各波段遥感数据进行几何校正。

2.2 ER-Mapper7.0遥感软件几何校正介绍

ER-Mapper7.0是一个大型专业影像处理软件, 由著名的澳大利亚ERM (Earth Resource Mapping) 公司研制, 软件发展了近20年, 遍及全球。几何校正功能位于软件主菜单“Process”中“Geocoding Wizard”项中, 或工具条“Ortho Geocoding Wizard”项中。其过程包括5大步骤: (1) 坐标变换方法选择, 如三角测量法、多项式法等; (2) 坐标变换参数选择, 如多项式法的n值取1, 2或3等; (3) 标准图件文件设置, 对应图像对图像、图像对地图等不同文件类型的几何校正。如图像对地图校正, 地图数据需转换为ER-Mapper7.0认可的tif、jpg等文件类型; (4) 编辑和建立控制点文件, 形成包括控制点坐标、精度等信息的数据文件, 文件后缀为gcp, 可以保存和多次使用; (5) 给定校正后图像文件名, 文件后缀为ers, 执行几何校正[7]。其中第 (4) 步是最为烦琐和至关重要的一步, 它关系到几何校正的成败和精度。

2.3 广东省从化市Landsat 7 ETM+几何校正

为使几何校正工作的顺利进行, 在开始前, 需要完成一项非常必要的准备工作, 就是反复视读ETM遥感图像和1:50万水系图件, 寻找和确定用于几何校正的控制点。由于广东省从化市ETM图像和1:50万水系图不是同时间采集的数据图件, 且比例尺也不一样, 这使从化市水系分布状况在这两个图件上会有许多不同点, 所以在选择控制点时, 需要遵循几何校正机理中谈到的控制点GCP采集要求。如果校正前, 熟悉并了解了不同图件上的水系特征信息, 做到可采集控制点的位置、分布和数量等都胸有成竹, 那么随后的几何校正工作也会事半功倍。

其次在地理信息系统MAPGIS环境下, 将1:50万水系的点、线矢量数据文件转换为tif后缀的图像文件格式。

在ER-Mapper7.0软件中进行几何校正的5个步骤的参数选择如下: (1) 选择“多项式”坐标变换方法; (2) 选择多项式n值为1; (3) 打开1:50万水系图件文件; (4) 在ETM+遥感图像和水系图之间编辑控制点。根据几何校正的机理分析, 选择多项式n值为1, 对应的坐标变换是二元一次线性变换。控制点数目选为3×6=18, 且控制点均匀分布在ETM+遥感图像和水系图中, 如图4所示。控制点精度参数RMS值小于1, 如表1所示; (5) 给定保存校正后图像文件名, 执行校正过程, 图5为几何校正后的ETM+图像。

3 几何校正效果分析

遥感图像的几何校正效果分析与评价, 一方面是从几何校正的机理参数上加以考虑, 分析和检查控制点的选择、分布、数量以及RMS值是否符合几何校正理论需要的条件, 从定量上把握几何校正的运算精度。另一方面是看几何校正的实际效果, 从视觉上, 查看不同图件上特征线的配准是否达到了预计的吻合程度。

第一方面的参数评价在ER-Mapper7.0等软件中比较容易实现, 观察校正中, 图像窗口中控制点分布是否均匀, 控制点数目是否大于3 (n+1) (n+2) , 控制点编辑表中的RMS值是否都小于1既可。第二方面的实际效果评价需要进行如下操作:在ER-Mapper7.0软件的算法窗口中建立两个Surface层, 分别装载ETM+图像数据和1:50万水系图数据, 系统自动将校正ETM+图像与1:50万水系图叠加显示, 如图6所示。从叠加图上, 可以看出整个ETM+图像位于1:50万水系图的右偏上方。水系图位于ETM+图像的底层, 调整ETM+图像的透明度, 可以观察水系图与ETM+图像中水系的重叠情况。追踪一些特征的水系线, 进行局部放大, 可以观察到几何校正后不同图件上水系分布重叠的实际效果。如局部放大图像左下角的溪流河部分, 如图7所示, 横穿北东向的溪流河在ETM+遥感图像以深色显示, 在水系图上以深色线条显示, 两种水系在叠加图上吻合得相当好。图右上方的溪流水库在两个图件上吻合也很好, 如图8所示。放大其它特征线部位, 都有类似的效果。说明实际的几何校正效果良好, 达到了几何校正的目的。

4 结论

遥感图像几何校正的理论性很强, 在实际应用中, 需要深刻理解其机理, 对控制点的选择、控制点分布和数量以及RMS值要有很好的把握, 才能作出较好的校正效果, 才能使校正后的遥感图像具有标准图件的几何精度和坐标系统。

通过多次实验性尝试, 发现下列误操作会带来一些错误的校正效果:

① 控制点数目偏少, 校正图像会产生局部扩张或发散问题。

② 控制点分布不均, 校正图像会在缺少控制点的部位变形。

③ 控制点RMS值过大, 校正图像会产生整体变形。

④ 控制点指示对象, 在不同图件上对应的不是同一物, 尽管上述3种情况都没有, 但校正图像会发生整体错位或变形。

前3种误操作可以在几何校正定量参数的选择时加以控制。后一种误操作较难察觉, 需要进行细致的特征线配准检查, 重新修正控制点坐标, 改善校正精度。

总之, 遥感图像的几何校正, 需要进行“理论”和“实际”两方面的效果评价, 方可得出较好的校正效果图。

摘要：几何校正是遥感图像处理的一个重要环节, 是削弱遥感图像与地面真实形态差异的重要手段。以广东省从化市ETM+遥感数据和GIS数据为例, 较深入地分析遥感图像几何校正原理并介绍ER-Mapper7.0遥感软件的几何校正应用, 指出提高遥感图像几何校正精度的有效途径和效果分析。

关键词：ETM+图像,GIS数据,几何校正,效果分析

参考文献