CN103297775B - 基于联合码率分配和控制的多光谱图像压缩系统与方法 - Google Patents

Abstract

本发明公开了一种基于联合码率分配和控制的多光谱图像压缩系统与方法，主要解决现有技术对多光谱图像压缩时码率分配不均导致重构质量不稳定及JPEG-LS码率控制不准确的问题。其实现过程是：对多光谱图像进行KL变换，得到变换后分量；将变换后分量进行数据重组，得到重组后分量；先将重组后分量的所有行数据进行第一次量化步长调整和量化编码，再对最后若干行的量化步长进行二次迭代调整，并根据量化步长进行量化编码，得到码流。本发明通过分量数据重组和统一量化步长调整，将码率分配和码率控制结合，同时解决了码率分配和码率控制的问题，码率控制精确且小于目标码率，各谱段重构质量均匀稳定，整体性能最优，适用于各种星上多光谱图像处理及传输。

Description

基于联合码率分配和控制的多光谱图像压缩系统与方法

技术领域

本发明涉及图像处理技术领域，更进一步涉及一种多光谱图像有损压缩系统与方法，可用于各种星上多光谱数据处理及传输。

背景技术

多光谱图像是由光谱成像仪对同一地物在数十个谱段上成像获取的同时包含空间信息和光谱信息的三维数据立方体，被广泛应用于资源勘探、目标识别、环境保护等方面。但是多光谱图像数据量非常庞大，需要采用有效的压缩技术以便图像的存储和传输。

多光谱图像具有较强的谱间相关性，有效利用这种相关性，可以大幅提高压缩性能，因此通常采用谱间变换+二维编码的方法对其进行高效的压缩。比如Du等人在论文“HyperspectralImageCompressionUsingJPEG2000andPrincipalComponentAnalysis”（IEEEGeoscienceandRemoteSensingLetters，4（2）：201-205）中提出利用KL（Karhunen–Loève）变换去除谱间相关性，再利用JPEG2000对变换分量进行压缩；张荣等人在论文“一种多光谱遥感图象的近无损压缩方法”（中国图象图形学报，3（10）：823-326）中采用KL变换结合预测树，对多光谱图像进行近无损压缩。

目前常用的谱间变换方法有DCT（DiscreteCosineTransform）变换、小波变换和KL变换等。其中KL变换是理论上最优的变换方法，可以完全去除数据间的相关性，被广泛的用于去除多光谱图像的谱间相关性。Hao等人在“MatrixFactorizationsforReversibleIntegerMapping”（IEEETransactionsonSignalProcessing，49（10）：2314-2324）中提出了KL变换的整数因式分解方法，对KL变换进行推广，实现了整数到整数的可逆变换，使其支持无损压缩。但是现有技术在有损压缩中均采用浮点KL变换，导致在低压缩比下有损压缩性能较低。

在二维编码方法中，静止图像无损/近无损压缩标准JPEG-LS(JointPhotographicExpertsGroupLosslessandnear_losslesscompressionofconituous_tonestillimage)，由于其在中高码率下具有较好的性能，而且具有较低的编码复杂度，被广泛应用于遥感图像压缩。但是JPEG-LS算法存在一个缺点，就是只能控制重构图像的失真，却不能控制压缩码率。而在遥感应用中，码率可控是一个非常重要的指标，如果压缩码率超过了存储空间和链路带宽，就会数据溢出，导致解码错误。因此需要设计JPEG-LS的码率控制方法来满足图像压缩对码率控制的要求。另外对于多光谱图像压缩而言，变换后各分量的能量分布不均，如果给各分量分配相同的码率，会导致某些分量质量较高，而某些分量质量较差，显然整体不是最优的。因此在各分量之间如何分配目标码率也是影响多光谱图像重构质量的一个关键问题。

码率分配和码率控制本质上是同一个问题，都是根据信号的能量信息来调整量化步长，使信号的总码率达到目标码率，并且整体质量达到最优。而现有技术通常将码率分配和码率控制两个问题分开处理，即首先根据能量分布将总目标码率分配给各个分量，然后在编码时分别调整每个分量的量化步长使其接近各分量的目标码率。这样做的缺点是，一方面由于首先分配码率，再对各个分量的码率分别控制，使各分量码率控制的误差会累积，导致整体码率控制不准确；另一方面在编码之前首先进行码率分配，获取各分量的先验知识较少，不能达到最优的码率分配。比如孙蕾等人在“基于谱间预测和码流预分配的高光谱图像压缩算法”（光学精密工程，2008，16(4):752-757）中提出根据各个分量的方差来分配相应的码率，这种方法虽然优于平均分配，但仍然导致各谱段的质量忽高忽低，不够稳定，个别谱段的重构质量太差，不能满足应用要求。

JPEG-LS的码率控制方法也得到了广泛的研究。西安空间无线电技术研究所在其专利申请“一种JPEG-LS图像压缩的码率控制方法”（专利申请号：201010617932.2，公开号：CN102088602A）以及西安电子科技大学在其硕士学位论文“JPEG-LS码率控制算法研究”中提出的码率控制方法均对图像进行分区，通过已编码区域的累积码率和目标码率的偏差量，对待编码区域的量化步长进行动态调整，使所有区域的码率之和达到总目标码率，但是码率控制的不够精确，而且实际码率在目标码率附近上下浮动，仍然会导致数据溢出，不能满足固定链路带宽的要求。

发明内容

本发明的目的在于提出一种基于联合码率分配和控制的多光谱图像压缩系统与方法，一是解决KL变换后各分量之间码率分配的问题，避免已有技术码率分配和码率控制分开，导致重构质量不稳定的缺点；二是解决对KL变换后分量进行JPEG-LS编码的码率控制问题，避免已有码率控制技术的实际码率在目标码率附近上下浮动，导致数据溢出的缺点。

实现本发明的技术关键在于KL变换后各分量的码率分配和码率控制，通过分析变换前后分量失真之间的关系，找到重构谱段失真相当的条件下各分量量化步长的比例关系；然后在码率可控JPEG-LS编码中，根据量化步长的关系将各个分量的数据重新组织在一起，对各分量的量化步长统一进行一次调整和二次迭代调整，使得码率分配和码率控制有效统一起来，并使得实际总码率逼近目标总码率。其技术方案如下：

一、本发明基于联合码率分配和控制的多光谱图像压缩系统，包括：

KL变换单元，用于对输入的N谱段多光谱图像{B¹,B²,…,B^k,…,B^N}进行KL变换，得到变换后分量{Y¹,Y²,…,Y^k,…,Y^N}，k＝1,2,…,N；

码率可控JPEG-LS编码单元，用于对变换后分量{Y¹,Y²,…,Y^k,…,Y^N}进行压缩编码，该单元包括：

数据重组模块，用于对变换后分量{Y¹,Y²,…,Y^k,…,Y^N}中的数据进行重新组织，得到重组后分量Y_CIL；

量化步长调整模块，用于对重组后分量Y_CIL的每一行数据进行量化步长动态调整，得到量化步长Near；

JPEG-LS编码模块，用于根据量化步长Near对重组后分量Y_CIL的每一行数据进行JPEG-LS编码，得到码流。

二、一种基于联合码率分配和控制的多光谱图像压缩方法，包括如下步骤：

（1）输入多光谱图像{B¹,B²,…,B^k,…,B^N}和目标压缩比Ratio，其中，输入的多光谱图像谱段数为N，宽为W，高为H，比特精度为I_p，Ratio＞0，k＝1,2,…,N；

（2）对输入多光谱图像{B¹,B²,…,B^k,…,B^N}进行KL变换，得到变换后分量{Y¹,Y²,…,Y^k,…,Y^N}；

（3）将变换后分量{Y¹,Y²,…,Y^k,…,Y^N}中的数据进行重新组织，得到重组后分量Y_CIL＝{Y_1,CIL,…,Y_i,CIL,…,Y_H,CIL}，其中每一行数据Y_i,CIL是按照分量序号k的顺序依次将各分量的第i行数据拼接构成，即i＝1,2,…,H；

（4）对重组后分量Y_CIL进行码率可控JPEG-LS编码：

（4a）参数设置及初始化：设置目标码率R_T=I_p/Ratio；初始化量化步长Near，若Ratio≤4，则Near=0；否则Near=4；设置码率门限a和Near调整步长b构成的正数对集合为：{(a⁺,b⁺)}={(0.5,10),(0.45,8),(0.35,6),(0.3,5),(0.25,3),(0.15,2),(0.05,1)}和负数对集合为：{(a^-,b^-)}＝{(-0.53,10),(-0.46,8),(-0.39,6),(-0.25,5),(-0.15,3),(-0.13,2),(-0.05,1)}；初始化迭代次数Times=0和迭代次数门限T_max=40，行索引i=1；

（4b）根据Near值对重组后分量Y_CIL中的第i行数据进行JPEG-LS量化编码，得到已编码数据的累积码率R_c；

（4c）计算已编码数据的累积码率R_c和目标码率R_T的差值slice_rate，通过依次比较slice_rate与码率门限值，对Near值进行一次调整，调整方法如下：

（4c1）判断slice_rate是否大于0，若slice_rate＞0，执行步骤（4c2）,否则执行步骤（4c3）；

（4c2）将slice_rate依次与正数对集合{(a⁺,b⁺)}={(0.5,10),(0.45,8),(0.35,6),(0.3,5),(0.25,3),(0.15,2),(0.05,1)}的码率门限a⁺比较，若slice_rate＞a⁺，则将Near值增加相应的b⁺，即将Near值调整为Near+b⁺；否则继续将slice_rate与下一个数对进行比较，直到最后一个数对；若slice_rate均小于所有正数对中的码率门限a⁺，则Near值保持不变；

（4c3）将slice_rate依次与负数对集合{(a^-,b^-)}＝{(-0.53,10),(-0.46,8),(-0.39,6),(-0.25,5),(-0.15,3),(-0.13,2),(-0.05,1)}中的码率门限比较，若slice_rate＜a^-，则将Near值减去相应的b^-，即将Near值调整为Near-b^-；否则继续将slice_rate与下一个数对进行比较，直到最后一个数对；若slice_rate均大于所有负数对中的码率门限a^-，则Near值保持不变；

（4d）判断第i行是否为倒数第2至倒数第P+1行中的某一行，如果是，存储第i行一次调整后的Near值；否则执行步骤（4e）；

（4e）判断第i行是否为倒数第P+1行，如果是，则存储第i行的上下文信息CX_P；否则执行步骤（4f）；

（4f）判断第i行是否为最后一行，如果不是，将i加1，执行步骤（4b）；如果是，执行步骤（4g）；

（4g）判断slice_rate是否大于0，如果是，执行步骤（4h），如果不是，则结束编码；

（4h）令i＝i-P，返回至倒数第P行，对最后P行的量化步长进行二次迭代调整，调整方法如下：

（4h1）将当前上下文信息CX更新为步骤（4e）存储的上下文信息CX_P，即CX=CX_P，并将Times加1；

（4h2）将第i行的Near值设置为步骤（4d）存储的第i行一次调整后的Near值，并将该Near值更新为Near+Times；

（4h3）根据该更新后的Near值对第i行数据进行JPEG-LS量化编码；

（4h4）判断第i行是否为最后一行，如果不是，将i加1，执行步骤（4h2）；如果是，判断Times是否大于迭代次数门限T_max，如果大于，则结束编码；否则，执行步骤（4g）。

（5）组织码流输出，编码结束。

本发明与现有技术相比较，具有如下优点：

第一，本发明将码率分配和码率控制有效统一起来，通过将变换后各分量重新组织在一起，统一对各分量进行量化步长调整和JPEG-LS编码，使得重构谱段的质量非常稳定，达到整体性能最优，避免了已有技术先分配各分量的目标码率，再分别利用调整量化步长控制各分量压缩码率所存在的重构质量不均衡的问题。

第二，本发明提出了基于二次迭代调整的JPEG-LS的动态精确码率控制技术，在完成所有行数据的第一次量化步长调整后，再对最后若干行数据的量化步长进行二次迭代调整，使得实际码率在不超过目标码率的前提下，十分逼近目标码率，克服了现有技术中码率控制不准确，上下浮动大的问题，保证了码率的精确控制，可以满足固定链路带宽的要求。

第三，本发明在低倍压缩时采用整数KL变换，在中高倍压缩时采用浮点KL变换，使得在各种压缩比时均具有较高的压缩性能；

附图说明

图1是本发明的压缩系统框图；

图2是本发明的压缩方法流程图；

图3是本发明KL变换前后分量量化失真关系示意图；

图4是本发明中基于二次迭代调整的码率可控JPEG-LS编码子流程图；

图5是本发明仿真实验的4组IKONOS多光谱图像；

图6是本发明的码率控制仿真结果曲线图；

图7是本发明的压缩性能仿真结果曲线图。

具体实施方式

参考图1，本发明基于联合码率分配和控制的多光谱图像压缩系统，包括KL变换单元和码率可控JPEG-LS编码单元。

所述的KL变换单元，对输入的多光谱图像{B¹,B²,…,B^k,…,B^N}进行KL变换，得到变换后分量{Y¹,Y²,…,Y^k,…,Y^N}，k＝1,2,…,N，输出给码率可控JPEG-LS编码单元，其中，输入的多光谱图像谱段数为N，宽为W，高为H，比特精度为I_p；

所述的码率可控JPEG-LS编码单元，包括数据重组模块、量化步长调整模块和JPEG-LS编码模块。该数据重组模块对变换后分量{Y¹,Y²,…,Y^k,…,Y^N}进行重新组织，得到重组后分量Y_CIL＝{Y_1,CIL,…,Y_i,CIL,…,Y_H,CIL}输出给量化步长调整模块和JPEG-LS编码模块，其中 为Y^k中的第i行数据，i＝1,2,…,H；该量化步长调整模块对重组后分量Y_CIL的每一行数据进行量化步长动态调整，得到量化步长Near输出给JPEG-LS编码模块；该JPEG-LS编码模块根据量化步长Near对重组后分量Y_CIL的每一行数据进行JPEG-LS编码，得到码流。

参照图2，本发明基于联合码率分配和控制的多光谱图像压缩方法详细步骤如下：

步骤1，输入多光谱图像{B₁,B₂,…,B_N}和目标压缩比Ratio，其中，输入的多光谱图像谱段数为N，宽为W，高为H，比特精度为I_p，Ratio＞0；

步骤2，对输入的多光谱图像{B¹,B²,…,B^k,…,B^N}进行KL变换，得到变换后分量{Y¹,Y²,…,Y^k,…,Y^N}，k＝1,2,…,N：

2.1）计算第k个谱段图像B^k中像素值的均值m_k：

$m_k=\left(\underset{i=1}{\overset{H}{\mathrm{\Σ}}}\underset{j=1}{\overset{W}{\mathrm{\Σ}}}{B}_{\mathrm{ij}}^k\right)/(W\×H),k=\mathrm{1,2},...,N,---1)$

其中，是第k个谱段图像中第i行第j列的像素值；

2.2）计算输入多光谱图像{B¹,B²,…,B^k,…,B^N}的谱间协方差矩阵C＝{C_u,v}，其中C_u,v为第u个谱段B^u和第v个谱段B^v的协方差，u＝1,2,…,N，v＝1,2,…,N：

C_u,v＝E{(B^u-m_u)(B^v-m_v)^T}＝E{B^uB^v}-m_um_v，2)

其中(·)^T表示矩阵的转置；E(·)表示数学期望，m_u和m_v分别是谱段B^u和B^v的均值；

2.3）对谱间协方差矩阵C进行正交相似变换，得到C的特征向量组{Φ₁,Φ₂,…,Φ_N}，用该特征向量组构成变换矩阵M＝[Φ₁Φ₂…Φ_N]^T，

所述的正交相似变换方法包括雅可比方法和QR分解方法等，本实施例中采用但不局限于隐式对称QR分解对C进行QR分解；

2.4）将输入的第k个谱段图像B^k按照行扫描的顺序存储为一维矢量B′^k，并对B′^k进行电平移位，得到移位后谱段X^k＝B′^k-m_k；

2.5）判断目标压缩比Ratio的数值，如果目标压缩比Ratio>2，执行步骤2.6），否则执行步骤2.7）；

2.6）对移位后谱段进行浮点KL变换，即将移位后谱段矩阵X＝{X¹,…,X^k,…X^N}与变换矩阵M相乘，得到变换后分量矩阵Y′＝{Y′¹,…,Y′^k,…,Y′^N}＝MX，执行步骤2.9）；

2.7）对变换矩阵M进行整数因式分解，得到转置矩阵P、第一下三角基本可逆矩阵L、上三角基本可逆矩阵U和第二下三角基本可逆矩阵S，即M＝PLUS；

2.8）对移位后谱段进行整数KL变换，即利用上述矩阵P、L、U、S通过多级提升分解，实现X→Y的整数到整数的可逆变换，得到变换后分量矩阵Y′＝{Y′¹,…,Y′^k,…,Y′^N}＝PLUSX；

2.9）将变换后矩阵Y′的第k行数据Y′^k按照行扫描的顺序存储为与输入多光谱图像大小相同的二维矩阵，作为变换后分量Y^k，该二维矩阵的宽为W，高为H。

步骤3，对变换后分量{Y¹,Y²,…,Y^k,…,Y^N}重新组织，得到重组后分量Y_CIL。

3a）重组原理分析

对于JPEG-LS量化编码而言，分配不同分量的目标码率相当于分配不同分量的量化步长。因此本发明通过分析各谱段重构质量相当条件下各分量量化步长的大小比例关系，将各分量{Y¹,Y²,…,Y^N}的数据重新组织，再按照其量化步长的比例关系，统一调整量化步长，达到同时进行码率分配和码率控制的目的。

根据KL变换公式 $\left(\begin{array}{c}{Y}^{\′1}\\ Y^{\′2}\\ \·\\ \·\\ \·\\ Y^{\′N}\end{array}\right)=M\left(\begin{array}{c}{B}^{\′1}-m_1\\ B^{\′2}-m_2\\ \·\\ \·\\ \·\\ B^{\′N}-m_N\end{array}\right),$ 对变换后分量矩阵{Y′¹,Y′²,…,Y′^N}直接进行反变换，得到原始谱段矩阵{B′¹,B′²,…,B′^N}如公式3）所示，其中M^T为变换矩阵M的转置：

$\left(\begin{array}{c}{B}^{\′1}\\ B^{\′2}\\ \·\\ \·\\ \·\\ B^{\′N}\end{array}\right)=M^T\left(\begin{array}{c}{Y}^{\′1}\\ Y^{\′2}\\ \·\\ \·\\ \·\\ Y^{\′N}\end{array}\right)+\left(\begin{array}{c}{m}_1\\ m_2\\ \·\\ \·\\ \·\\ m_N\end{array}\right)---3)$

对分量{Y′¹,Y′²,…,Y′^N}进行量化、编码及重构得到重构分量再对重构分量进行反变换得到重构谱段如公式4）所示：

$\left(\begin{array}{c}{\tilde{B}}^{\′1}\\ {\tilde{B}}^{\′2}\\ \·\\ \·\\ \·\\ {\tilde{B}}^{\′N}\end{array}\right)=M^T\left(\begin{array}{c}{\tilde{Y}}^{\′1}\\ {\tilde{Y}}^{\′2}\\ \·\\ \·\\ \·\\ {\tilde{Y}}^{\′N}\end{array}\right)+\left(\begin{array}{c}{m}_1\\ m_2\\ \·\\ \·\\ \·\\ m_N\end{array}\right)---4)$

将公式3）和公式4）直接相减，即可得到重构谱段失真ΔB＝{ΔB¹,…,ΔB^k,…,ΔB^N}和分量失真ΔY＝{ΔY¹,…,ΔY^k,…,ΔY^N}的关系如公式5）所示,其中ΔB^k是第k个谱段B′^k的失真，ΔY^k是第k个分量Y′^k的失真：

$\left(\begin{array}{c}{\mathrm{\ΔB}}^1\\ {\mathrm{\ΔB}}^2\\ \·\\ \·\\ \·\\ {\mathrm{\ΔB}}^N\end{array}\right)=M^T\left(\begin{array}{c}{\mathrm{\ΔY}}^1\\ {\mathrm{\ΔY}}^2\\ \·\\ \·\\ \·\\ {\mathrm{\ΔY}}^N\end{array}\right)---5)$

根据公式5）分析在各重构谱段失真{ΔB^k}相当时，各分量的量化步长{Q^k}的比例关系。在JPEG-LS中，量化失真是由量化步长决定的。如果量化步长为Q，那么该像素的失真|Δ|≤Q。以二维向量{Y′¹,Y′²}为例，如图3所示，若Y′¹、Y′²的量化步长都为Q，则分量的失真ΔY¹和ΔY²分布在图3中的正方形内。M^T是正交矩阵，根据正交矩阵的性质，正交矩阵与向量相乘，只改变向量的角度，而不改变向量的长度，相当于对向量进行旋转，即M^T作用于分量失真ΔY，相当于对该正方形进行旋转，得到谱段失真ΔB，如公式6）所示

||ΔB||＝||M^TΔY||＝||ΔY||.6）

经过旋转后，对于大多数像素而言，即半径为Q的圆形阴影区域内的像素，ΔB仍然在正方形内，即B^k的失真仍然限制在量化步长Q内；而对其余少量像素而言，即半径为的圆内像素，其失真也限制在以内，因此每个谱段的失真都是相当的。对于更多维数的情况也是类似的。

通过上述分析可以得到以下结论：当各分量采用相同的量化步长Q^k＝Q时，各重构谱段的失真也是相当的。

因此，需要对变换后分量进行重新组织，将各分量{Y¹,Y²,…,Y^k,…,Y^N}的数据组织在一起，得到重组后分量Y_CIL，使得在后续码率可控JPEG-LS编码中可以对Y_CIL中各分量数据统一调整量化步长，这样构成Y_CIL分量的各分量{Y¹,Y²,…,Y^k,…,Y^N}的调整量化步长始终是相等的，可以保证各谱段的重构质量相当。

3b）将变换后分量{Y¹,Y²,…,Y^k,…,Y^N}按照CIL（Component-Interleaved-by-Line）的格式拼接在一起，即类似于多光谱图像的BIL（Band-Interleaved-by-Line）数据组织方式进行重新组织，重组步骤如下：

3b1）设置行索引i=1，分量索引k=1；

3b2）将第k个分量的第i行数据依次存储作为重组后分量Y_CIL的第i行数据Y_i,CIL；

3b3）令i=i+1，若i≤H,执行步骤3b2），否则令i=1，执行步骤3b4）；

3b4）令k＝k+1，若k≤N，执行步骤3b2）；否则重组结束，得到重组后分量Y_CIL＝{Y_1,CIL,…,Y_i,CIL,…,Y_H,CIL}，其中Y_CIL的大小为NW×H，其中，N为输入多光谱图像的谱段数，W为多光谱图像的宽，H为多光谱图像的高。

步骤4，对重组后分量Y_CIL进行基于二次迭代调整的码率可控JPEG-LS编码，得到压缩码流。

参照图4，本步骤具体实现如下：

（4.1）参数设置及初始化：设置目标码率R_T=I_p/Ratio；初始化量化步长Near，设置码率门限和Near调整步长门限。在本实施例中，若Ratio≤4，则Near=0；否则Near=4；设置码率门限a和Near调整步长b构成的正数对集合为：{(a⁺,b⁺)}={(0.5,10),(0.45,8),(0.35,6),(0.3,5),(0.25,3),(0.15,2),(0.05,1)}，负数对集合为：{(a^-,b^-)}＝{(-0.53,10),(-0.46,8),(-0.39,6),(-0.25,5),(-0.15,3),(-0.13,2),(-0.05,1)}；初始化迭代次数Times=0和迭代次数门限T_max=40，行索引i=1。

（4.2）根据Near值对重组后分量Y_CIL中的第i行数据进行JPEG-LS量化编码，得到已编码数据的累积码率R_c；

（4.3）计算已编码数据的累积码率R_c和目标码率R_T的差值slice_rate，通过依次比较slice_rate与码率门限值，对Near值进行一次调整，调整方法如下：

（4.3a）判断slice_rate是否大于0，若slice_rate＞0，执行步骤（4.3b）,否则执行步骤（4.3c）；

（4.3b）将slice_rate依次与正数对集合{(a⁺,b⁺)}={(0.5,10),(0.45,8),(0.35,6),(0.3,5),(0.25,3),(0.15,2),(0.05,1)}的码率门限a⁺比较，若slice_rate＞a⁺，则将Near值增加相应的b⁺，即将Near值调整为Near+b⁺，否则继续将slice_rate与正数对集合的下一个数对进行比较，直到最后一个数对；若slice_rate均小于所有正数对中的码率门限a⁺，则Near值保持不变；

（4.3c）将slice_rate依次与负数对集合{(a^-,b^-)}＝{(-0.53,10),(-0.46,8),(-0.39,6),(-0.25,5),(-0.15,3),(-0.13,2),(-0.05,1)}中的码率门限a^-比较，若slice_rate＜a^-，则将Near值减去相应的b^-，即将Near值调整为Near-b^-；否则继续将slice_rate与负数对集合中的下一个数对进行比较，直到最后一个数对；若slice_rate均大于所有负数对中的码率门限a^-，则Near值保持不变；

（4.4）判断第i行是否为倒数第2至倒数第P+1行中的某一行，如果是，存储第i行一次调整后的Near值；否则执行步骤（4.5），在本实施例中P=4；

（4.5）判断第i行是否为倒数第P+1行，如果是，则存储第i行的上下文信息CX_P，否则执行步骤（4.6）；

（4.6）判断第i行是否为最后一行，如果不是，将i加1，执行步骤（4.2）；如果是，执行步骤（4.7）；

（4.7）判断slice_rate是否大于0，如果是，执行步骤（4.8），如果不是，则结束编码，得到压缩码流；

（4.8）令i＝i-P，返回至倒数第P行，对最后P行的量化步长进行二次迭代调整：

（4.8a）将当前上下文信息CX更新为步骤（4.5）存储的上下文信息CX_P，即CX=CX_P，并将Times加1；

（4.8b）将第i行的Near值设置为步骤（4.4）存储的第i行一次调整后Near值，并将该Near值更新为Near+Times；

（4.8c）根据该更新后的Near值对当前第i行数据进行JPEG-LS量化编码；

（4.8d）判断第i行是否为最后一行，如果不是，将i加1，执行步骤（4.8b）；如果是，执行步骤（4.8e）；

（4.8e）判断Times是否大于迭代次数门限T_max，如果大于，则结束编码，得到压缩码流；否则，返回步骤（4.7）。

步骤5，组织码流并输出，编码结束。

（5.1）将编码必要参数组织为码流头信息并输出，所述编码必要参数包括：每个谱段的像素均值m_k、KL变换矩阵M、最后P行二次迭代调整后的Near值；

（5.2）输出重组后分量Y_CIL的压缩码流。

本发明的效果可以通过以下仿真实验做进一步说明。

（1）仿真条件

软件采用MicrosoftVisualC++6.0集成开发软件和C语言在Microsoft公司的Windows7环境上实现。本发明选取了4组IKONOS多光谱图像，每组多光谱图像包含红、绿、蓝和近红外4个谱段，图像的空间大小为1024×1024像素，每个像素采用10比特存储，如图5所示，其中图5（a）表示第1组多光谱图像的第一个谱段、图5（b）表示第2组多光谱图像的第一个谱段、图5（c）表示第3组多光谱图像的第一个谱段、图5（d）表示第4组多光谱图像的第一个谱段。

（2）仿真内容

（2.1）码率控制性能仿真

选用了两种方法对本发明提出的基于二次迭代调整的JPEG-LS码率控制方法的性能进行评估。这两种方法分别是西安电子科技大学硕士学位论文“JPEG-LS码率控制算法研究”中提出的基于分区一次调整的码率控制方法和本发明方法提出的基于二次迭代调整的JPEG-LS码率控制方法，仿真结果如图6所示，图6中的横坐标1-1、1-2、…、4-4分别对应了上述4组多光谱图像的4个谱段，其中：

图6（a）为用上述两种方法以3倍独立压缩多光谱图像各谱段时的实际码率；图6（b）为用上述两种方法以3倍独立压缩多光谱图像各谱段时的峰值信噪比PSNR。

从图6（a）中可以看出，本发明与分区一次调整方法相比，控制码率更加精确，本发明的实际码率在低于目标码率的前提下，非常逼近目标码率，即与目标码率的偏差仅仅为0.0064。

从图6（b）可以看出，由于本发明码率控制更加准确，重构质量也明显优于分区一次调整，平均PSNR比分区一次调整提高了大约1dB。

表1显示了上述两种方法的编码运行时间，从表1中可以看出本发明提出的码率控制方法仅仅对少数行进行了简单的二次迭代调整，因此并没有增加太多编码复杂度，编码运行时间仅仅比分区一次调整方法增加了大约10ms。

表1：编码时间比较

	分区一次调整	本发明
			编码时间（ms）	510.38	521.63

（2.2）压缩性能仿真

（2.2a）仿真方法

仿真方法包含本发明和两种现有方法，其中第一种现有方法是基于码率可控JPEG-LS的独立压缩方法，第二种现有方法是基于KL变换和码率可控JPEG-LS的联合压缩方法，并采用孙蕾等人在“基于谱间预测和码流预分配的高光谱图像压缩算法”中提出的基于方差的码率分配方法。

（2.2b）仿真内容

利用上述三种方法在2倍和3倍压缩比下对图5中的4组多光谱图像进行压缩，并用每组多光谱图像的PSNR平均值和PSNR标准差作为压缩性能的评价指标，为了便于统计平均值和标准差，将无损时的PSNR记为100dB，结果如图7。其中：

图7（a）为用上述3种方法在2倍压缩下各组多光谱图像的PSNR平均值；

图7（b）为用上述3种方法在3倍压缩下各组多光谱图像的PSNR平均值；

图7（c）为用上述3种方法在2倍压缩下各组多光谱图像的PSNR标准差；

图7（d）为用上述3种方法在3倍压缩下各组多光谱图像的PSNR标准差。

图7中的横坐标1、2、3、4对应了上述4组多光谱图像。

从图7（a）和图7（b）可以看出，本发明的重构质量PSNR平均值明显高于独立压缩方法和基于方差的码率分配方法，整体压缩性能最优，尤其是图7（a）所示的2倍压缩，本发明有3组多光谱图像达到了无损压缩，即100dB。

从图7（c）和图7（d）可以看出，本发明系统重构的多光谱图像PSNR标准差明显小于独立压缩方法和基于方差的码率分配方法，本发明恢复的各谱段PSNR的标准差均保持在0.1以内，因此重构质量的稳定性最好，避免了其它方法重构谱段质量忽高忽低，不够稳定的问题。

以上描述仅是本发明的一个具体实例，不构成对本发明的任何限制。显然对于本领域的专业人员来说，在了解了本发明内容和原理后，都可能在不背离本发明原理、结构的情况下，进行形式和细节上的各种修正和改变，但是这些基于本发明思想的修正和改变仍在本发明的权利要求保护范围之内。

Claims (2)

1.一种基于联合码率分配和控制的多光谱图像压缩方法，包括如下步骤：

(1)输入多光谱图像{B¹,B²,…,B^k,…,B^N}和目标压缩比Ratio，其中，输入的多光谱图像谱段数为N，宽为W，高为H，比特精度为I_p，Ratio>0，k＝1,2,…,N；

(2)对输入多光谱图像{B¹,B²,…,B^k,…,B^N}进行KL变换，得到变换后分量{Y¹,Y²,…,Y^k,…,Y^N}；

(3)将变换后分量{Y¹,Y²,…,Y^k,…,Y^N}中的数据进行重新组织，得到重组后分量Y_CIL＝{Y_1,CIL,…,Y_i,CIL,…,Y_H,CIL}，其中每一行数据Y_i,CIL是按照分量序号k的顺序依次将各分量的第i行数据{Y_i ^k}拼接构成，即Y_i,CIL＝{Y_i ¹Y_i ²…Y_i ^k…Y_i ^N}，i＝1,2,…,H；

(4)对重组后分量Y_CIL进行码率可控JPEG-LS编码：

(4a)参数设置及初始化：设置目标码率R_T＝I_p/Ratio；初始化量化步长Near，若Ratio≤4，则Near＝0；否则Near＝4；设置码率门限a和Near调整步长b构成的正数对集合为：{(a⁺,b⁺)}＝{(0.5,10),(0.45,8),(0.35,6),(0.3,5),(0.25,3),(0.15,2),(0.05,1)}和负数对集合为：{(a^-,b^-)}＝{(-0.53,10),(-0.46,8),(-0.39,6),(-0.25,5),(-0.15,3),(-0.13,2),(-0.05,1)}；初始化迭代次数Times＝0和迭代次数门限T_max＝40，行索引i＝1；

(4b)根据Near值对重组后分量Y_CIL中的第i行数据进行JPEG-LS量化编码，得到已编码数据的累积码率R_c；

(4c)计算已编码数据的累积码率R_c和目标码率R_T的差值slice_rate，通过依次比较slice_rate与码率门限值，对Near值进行一次调整，调整方法如下：

(4c1)判断slice_rate是否大于0，若slice_rate>0，执行步骤(4c2),否则执行步骤(4c3)；

(4c2)将slice_rate依次与正数对集合{(a⁺,b⁺)}＝{(0.5,10),(0.45,8),(0.35,6),(0.3,5),(0.25,3),(0.15,2),(0.05,1)}的码率门限a⁺比较，若slice_rate>a⁺，则将Near值增加相应的b⁺，即将Near值调整为Near+b⁺；否则继续将slice_rate与下一个数对进行比较，直到最后一个数对；若slice_rate均小于所有正数对中的码率门限a⁺，则Near值保持不变；

(4c3)将slice_rate依次与负数对集合{(a^-,b^-)}＝{(-0.53,10),(-0.46,8),(-0.39,6),(-0.25,5),(-0.15,3),(-0.13,2),(-0.05,1)}中的码率门限比较，若slice_rate<a^-，则将Near值减去相应的b^-，即将Near值调整为Near-b^-；否则继续将slice_rate与下一个数对进行比较，直到最后一个数对；若slice_rate均大于所有负数对中的码率门限a^-，则Near值保持不变；

(4d)判断第i行是否为倒数第2至倒数第P+1行中的某一行，如果是，存储第i行一次调整后的Near值；否则执行步骤(4e)；

(4e)判断第i行是否为倒数第P+1行，如果是，则存储第i行的上下文信息CX_P；否则执行步骤(4f)；

(4f)判断第i行是否为最后一行，如果不是，将i加1，执行步骤(4b)；如果是，执行步骤(4g)；

(4g)判断slice_rate是否大于0，如果是，执行步骤(4h)，如果不是，则结束编码；

(4h)令i＝i-P，返回至倒数第P行，对最后P行的量化步长进行二次迭代调整，调整方法如下：

(4h1)将当前上下文信息CX更新为步骤(4e)存储的上下文信息CX_P，即CX＝CX_P，并将Times加1；

(4h2)将第i行的Near值设置为步骤(4d)存储的第i行一次调整后的Near值，并将该Near值更新为Near+Times；

(4h3)根据该更新后的Near值对第i行数据进行JPEG-LS量化编码；

(4h4)判断第i行是否为最后一行，如果不是，将i加1，执行步骤(4h2)；如果是，判断Times是否大于迭代次数门限T_max，如果大于，则结束编码；否则，执行步骤(4g)；

(5)组织码流输出，编码结束。

2.根据权利要求1所述基于联合码率分配和控制的多光谱图像压缩方法，其特征在于步骤(2)所述的对输入的多光谱图像{B¹,B²,…,B^k,…,B^N}，k＝1,2,…N，进行KL变换，按如下步骤进行：

(2a)计算第k个谱段图像B^k中像素值的均值m_k；

(2b)计算多光谱图像{B¹,B²,…,B^k,…,B^N}的谱间协方差矩阵C；

(2c)对谱间协方差矩阵C进行正交相似变换，得到C的特征向量{Φ₁,Φ₂,…,Φ_k,…,Φ_N}，构成变换矩阵M＝[Φ₁Φ₂…Φ_k…Φ_N]^T；

(2d)将输入的第k个谱段图像B^k按照行扫描的顺序存储为一维矢量B′^k，对B′^k进行电平移位，得到移位后谱段X^k＝B′^k-m_k；

(2e)判断目标压缩比Ratio的数值，如果目标压缩比Ratio>2，则执行步骤(2f)，否则执行步骤(2g)；

(2f)利用变换矩阵M对移位后谱段进行浮点KL变换，即将移位后谱段矩阵X＝{X¹,…,X^k,…,X^N}与变换矩阵M相乘，得到变换后分量矩阵Y′＝{Y′¹,…,Y′^k,…,Y′^N}＝MX，执行步骤(2i)；

(2g)对变换矩阵M进行整数因式分解，将M分解为转置矩阵P、第一下三角基本可逆矩阵L、上三角基本可逆矩阵U和第二下三角基本可逆矩阵S，即M＝PLUS；

(2h)利用所述矩阵P、L、U、S对移位后谱段矩阵X进行整数KL变换，得到变换后分量矩阵Y′＝{Y′¹,…,Y′^k,…,Y′^N}＝PLUSX；

(2i)将变换后分量矩阵Y′的第k行数据Y′^k按照行扫描的顺序存储为与输入多光谱图像大小相同的二维矩阵作为变换后分量Y^k，即宽为W，高为H，变换结束。