CN109191404A - 一种基于e-3dtv正则的高光谱图像修复方法 - Google Patents
一种基于e-3dtv正则的高光谱图像修复方法 Download PDFInfo
- Publication number
- CN109191404A CN109191404A CN201811046032.XA CN201811046032A CN109191404A CN 109191404 A CN109191404 A CN 109191404A CN 201811046032 A CN201811046032 A CN 201811046032A CN 109191404 A CN109191404 A CN 109191404A
- Authority
- CN
- China
- Prior art keywords
- matrix
- 3dtv
- canonical
- formula
- follows
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 238000001228 spectrum Methods 0.000 title claims abstract description 93
- 238000000034 method Methods 0.000 title claims abstract description 75
- 239000011159 matrix material Substances 0.000 claims abstract description 87
- 230000003595 spectral effect Effects 0.000 claims abstract description 35
- 230000006835 compression Effects 0.000 claims abstract description 29
- 238000007906 compression Methods 0.000 claims abstract description 29
- 230000008439 repair process Effects 0.000 claims abstract description 13
- 230000000694 effects Effects 0.000 claims description 39
- 238000010586 diagram Methods 0.000 claims description 18
- 239000000284 extract Substances 0.000 claims description 13
- 238000010276 construction Methods 0.000 claims description 8
- 238000000354 decomposition reaction Methods 0.000 claims description 8
- 230000003190 augmentative effect Effects 0.000 claims description 6
- 238000004364 calculation method Methods 0.000 claims description 6
- 230000006870 function Effects 0.000 claims description 6
- 230000000630 rising effect Effects 0.000 claims description 6
- 238000005070 sampling Methods 0.000 claims description 6
- 238000004422 calculation algorithm Methods 0.000 claims description 5
- 230000008569 process Effects 0.000 claims description 5
- 238000009795 derivation Methods 0.000 claims description 4
- 238000012545 processing Methods 0.000 claims description 4
- 239000000047 product Substances 0.000 claims description 4
- 230000008859 change Effects 0.000 claims description 3
- 239000012467 final product Substances 0.000 claims description 3
- 238000005457 optimization Methods 0.000 claims description 3
- 230000009466 transformation Effects 0.000 claims description 3
- 230000008602 contraction Effects 0.000 claims 3
- 230000001010 compromised effect Effects 0.000 claims 2
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 claims 1
- 230000002708 enhancing effect Effects 0.000 abstract description 4
- 230000007547 defect Effects 0.000 abstract description 2
- 238000002474 experimental method Methods 0.000 description 5
- 235000008331 Pinus X rigitaeda Nutrition 0.000 description 4
- 235000011613 Pinus brutia Nutrition 0.000 description 4
- 241000018646 Pinus brutia Species 0.000 description 4
- 230000005540 biological transmission Effects 0.000 description 4
- 235000013399 edible fruits Nutrition 0.000 description 4
- 238000009499 grossing Methods 0.000 description 4
- 238000009412 basement excavation Methods 0.000 description 2
- 230000006872 improvement Effects 0.000 description 2
- 239000000463 material Substances 0.000 description 2
- 239000000126 substance Substances 0.000 description 2
- 230000003321 amplification Effects 0.000 description 1
- 239000005427 atmospheric aerosol Substances 0.000 description 1
- 230000008901 benefit Effects 0.000 description 1
- 230000019771 cognition Effects 0.000 description 1
- 230000008878 coupling Effects 0.000 description 1
- 238000010168 coupling process Methods 0.000 description 1
- 238000005859 coupling reaction Methods 0.000 description 1
- 238000013461 design Methods 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 238000000605 extraction Methods 0.000 description 1
- 230000004927 fusion Effects 0.000 description 1
- 238000003384 imaging method Methods 0.000 description 1
- 238000002156 mixing Methods 0.000 description 1
- 239000000203 mixture Substances 0.000 description 1
- 238000003199 nucleic acid amplification method Methods 0.000 description 1
- 230000003287 optical effect Effects 0.000 description 1
- VMXUWOKSQNHOCA-UKTHLTGXSA-N ranitidine Chemical compound [O-][N+](=O)\C=C(/NC)NCCSCC1=CC=C(CN(C)C)O1 VMXUWOKSQNHOCA-UKTHLTGXSA-N 0.000 description 1
- 238000011084 recovery Methods 0.000 description 1
- 238000012360 testing method Methods 0.000 description 1
- 230000007704 transition Effects 0.000 description 1
- 230000000007 visual effect Effects 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/70—Denoising; Smoothing
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/77—Retouching; Inpainting; Scratch removal
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10032—Satellite or aerial image; Remote sensing
- G06T2207/10036—Multispectral image; Hyperspectral image
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/20—Special algorithmic details
- G06T2207/20212—Image combination
- G06T2207/20224—Image subtraction
Landscapes
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Engineering & Computer Science (AREA)
- Theoretical Computer Science (AREA)
- Image Processing (AREA)
Abstract
一种基于E‑3DTV正则的高光谱图像修复方法，将含噪音的原三维高光谱数沿光谱维度展开成矩阵，初始化噪音项和待修复高光谱数据的矩阵表示以及ADMM框架下的其他模型变量及参数；对待修复的高光谱数据沿水平、竖直、光谱三个维度做差分操作，得到三个不同方向的梯度图，并沿光谱维度展开成矩阵；对三个方向的梯度图矩阵分别做低秩UV分解，对梯度图的基矩阵加稀疏性约束，得到E‑3DTV正则；对待修复的数据加E‑3DTV正则，写出优化模型，利用ADMM框架进行迭代求解，到迭代稳定，获得修复图像和噪音。本发明对高光谱图像数据进行去噪和压缩重构，通过对传统3DTV的增强改进，可兼顾梯度图像的结构相关性和稀疏性，克服传统的3DTV只能刻画梯度图的稀疏性，而忽略相关性的缺陷。
Description
技术领域
本发明涉及对传统的3DTV的一种增强改进，使得新的正则可以兼顾高光谱图梯度图像的结构的相关性和稀疏性，而不同于传统的3DTV只能刻画梯度图的稀疏性，而忽略了其本质的相关性。具体涉及针对高光谱图像数据的去噪和压缩重构方法。
背景技术
卫星上携带的光谱成像仪可以在同一时刻对同一地物形成许多连续谱段下的图像，因此在光谱维度每一个像素都会形成一条独特的光谱曲线，这是区分不同物质的一个显著特征。因此，高光谱图像相比传统图像在光谱维度具有更为丰富的信息，因此被广泛应用到地物识别，城市规划，军事勘测等应用问题中。
然而，地端到卫星之间的成像距离通常很长，成像中往往会受到大气气溶胶、光照、水气等因素的污染，使得高光谱图像中往往混杂大量噪音，从而使得部分高光谱图像的结构纹理模糊，物质的光谱曲线被严重破坏，这些现象都严重影响了针对高光谱图像进行应用与处理的后续目标，因此针对高光谱图像的修复具有非常重要的实用价值。除此之外，高光谱图像数据需要被传输到地面来进行进一步的分析与处理任务，由于高光谱图像通常蕴含高谱段的丰富信息，其数据量往往非常巨大。如何有效设计传输重构算法来快速有效的解决高光谱图像传输和储存问题，是针对高光谱图像另一个具有重要意义而挑战性的任务。由于近年来压缩重构理论的发展，压缩重构逐渐倾向于突破奈奎斯特采样的限制，所需传输的采样子个数期待大大减少，高光谱的传输问题也转变为利用压缩重构获得好的重构效果的问题。
以上高光谱的去噪和压缩重构问题为最具代表性的两类高光谱图像修复问题。
高光谱图像修复技术主要基于构建高光谱图像先验信息来进行恢复。根据其先验只是的分类，可以将高光谱图像修复方法分为三类。基于高光谱图像光谱之间的相关性先验知识构建的方法，基于高光谱图像的局部平滑性先验构造的方法，与基于高光谱图非局部自相似性先验构造的方法。实际方法可能会综合采用以上两类及多类先验来构造高光谱图像修复算法。
其中基于光谱相关性先验知识构造的方法，是指利用高光谱图像在光谱维度上的相关性，通过构造一组基及在基下的表达系数来表示高光谱图像。比较常见的方法包括低秩类方法，包括低秩矩阵和低秩张量类方法，以及字典表达类方法。此类方法均为基于高光谱的主体信息可以通过稀疏表达而噪音不能被稀疏表达的基本礼节来实现噪音和干净数据分离的效果。
基于图像的非局部自相似先验类的方法，是在高光谱图像块层面上进行的。将图片分成许多小块后，有许多相似的小块出现。这样就可以将图片的小块分成许多相似的小组，由于每个小组的小块都比较相似，可利用低秩和字典学习等方法对每组块进行共有信息的抽取，从而实现对每组块共同的修复，最终将所有修复组进行合并就可以最终得到修复后的高光谱图。
基于图像的局部平滑性先验类的方法，是基于如下统计规律而构建：相邻的像素一般是同类的物质，因此相邻像素的灰度值差异很小，在不同物质的连接处，像素的灰度值的差异很大，这种灰度值的差异反映了图像的主体结构，即是图像的边缘结构。可以用图片的梯度图的稀疏性来刻画这种灰度值差异的先验，因此一般用全变差正则(TV)来刻画这种先验信息。当对图像使用梯度或者差分算子操作时，梯度图可以显著的提取图像边缘结构，这对高光谱图像的修复会起到很好的作用。
现有的技术尽管是基于对高光谱图的先验知识来构造模型，但其先验并未能充分精细刻画高光谱图像的信息。如基于高光谱图的相似性的方法忽略了高光谱图的局部平滑性，而基于图像的局部平滑性的又不能很好的刻画光谱的相似性，尽管也有一些基于两类甚至三类先验融合的方法，但这些方法的先验知识通常采用简单的叠加，并不能充分的刻画高光谱图丰富的先验结构耦合的信息。因此寻找能够充分表达高光谱图信息先验及其先验编码方式对于推动高光谱图修复方法的性能具有非常重要的作用与意义。
发明内容
本发明的目的在于提供一种基于增强三维全变差正则(Enhanced 3-dimensionaltotal variation,简称E-3DTV)的高光谱图像(Hyperspectral image)修复方法，对高光谱图像数据进行去噪和压缩重构，通过对传统3DTV的增强改进，使得新的正则可以兼顾梯度图像的结构相关性和稀疏性，克服传统的3DTV只能刻画梯度图的稀疏性，而忽略了相关性的缺陷。
为了实现上述目的，本发明采用的技术方案是：
一种基于E-3DTV正则的高光谱图像修复方法，包括如下步骤：
步骤S1：将含噪音的原三维高光谱数据沿着光谱维度展开成矩阵得到Y∈Rmn×T，其中m,n代表每个谱段下图像的长和宽，T代表高光谱的光谱数，初始化噪音项ε∈Rm×n×T的矩阵形式E∈Rmn×T和待修复的高光谱数据的矩阵表示X∈Rmn×T以及ADMM框架下的其他模型变量及参数；
步骤S2：对待修复的高光谱数据沿着水平、竖直、光谱三个维度做差分操作，差分算子记为得到的三个不同方向的梯度图，并沿着光谱维度展开成矩阵，记水平、竖直、光谱方向的梯度图分别为i＝1,2,3；
步骤S3：对三个方向的梯度图矩阵分别做低秩UV分解，其中基矩阵U反应反映了梯度图的共有的结构和特征信息，V是系数矩阵，它反映了梯度图矩阵的每一列由基矩阵表达时的系数；对梯度图的基矩阵U加稀疏性约束，得到E-3DTV正则；
步骤S4：根据步骤S3，对待修复的X加E-3DTV正则，写出优化模型，利用ADMM框架进行迭代求解，到迭代稳定，获得修复图像X和噪音E。
所述步骤S1中，ADMM框架下的其他模型变量及参数包括ADMM中需要用到的上升乘子矩阵M和上升乘子常数的初始值μ和每次上升的倍率ρ，M＝0，μ的经验取值在[1e-3,1e-2]之间选择，ρ的经验取值在[1.05,1.8]之间选择。
所述步骤S2中三个方向的差分图计算公式定义如下所示：
其中，i,j,k分别为高光谱图像对应的沿空间长度、宽度与谱段方向的坐标位置，表示原高光谱图像在位置(i,j,k)的灰度值，s＝1,2,3分别表示的是原高光谱的在位置的灰度值表示沿着空间长度、宽度与谱段对高光谱图进行差分运算之后得到的梯度图在位置(i,j,k)的取值。分别表示对应高光谱图像的(i,j,k)位置沿着空间长度、宽度与谱段方向分别平移一个位置的像素灰度值。为了保证差分图和待修复图尺寸一致便于计算，对差分图做zero-padding处理，将梯度图沿着光谱维度展开为矩阵，得到如下的公式：
其中unfold表示沿着光谱方向展开。
所述步骤S3中，E-3DTV正则的形式如下：
其中，Ui为梯度矩阵的基，由于Ui中包含了梯度矩阵的共有结构或者特征，因此给Ui加上了一范数‖Ui‖1来约束基的稀疏性；Vi为梯度矩阵在基矩阵下的系数矩阵。为了保证上述优化问题解的唯一性，给系数矩阵加上ViVi T＝I的限制，这里I表示单位矩阵；r为低秩UV分解中的秩，r<min(mn,T)，蕴含了重建矩阵的低秩性信息。
所述步骤S4中，建立的优化模型如下所示：
s.t.Y＝X+E (4)
其中τ为调节参数，作用在于协调E-3DTV正则和误差项的作用折中，Y为输入的带噪音数据，E为噪音项，X为待修复的高光谱图数据，等式(4)展开即可得到下面的公式(5)：
所述步骤S4中，利用ADMM框架进行迭代求解的步骤如下：
S6.1)给出公式(5)的增广拉格朗日函数：
其中，Mi,Γ是拉格朗日乘子，μ是上升乘子常数，<.，.>表示内积运算，表示矩阵的F范数；
S6.2)建立交替方向乘子法的迭代格式与终止条件：
μk+1＝ρμk (12)
迭代终止条件为：
S6.3)对问题(6)、(7)、(8)、(9)进行求解，给出迭代的具体算式；
S6.4)设置迭代的初始值为：E0＝0，由奇异值分解方法作用在Y上产生；
S6.5)进行公式(6)-(12)的迭代运算，直到迭代满足终止条件。
通过所述步骤S1-S4实现高光谱图像去噪；在所述步骤S1-S4的基础上，增加一步由少量采样数据y来计算重构数据的过程，实现高光谱图像压缩重构。
所述高光谱图像压缩重构的具体流程如下：
S8.1)压缩重构问题简述
构造采样算子Ψ＝D.H.P，其中D是随机下采样算子，H是随机置换矩阵，P是WalshHadamard变换算子，固定采样算子Ψ之后，高光谱数据Z∈Rmn×T经Ψ的作用会得到采样过后的观测子y∈Rl,l＜＜mnT，压缩重构的目的就是通过观测y重构出Z。压缩重构蕴含内在等式约束为y＝Ψ(Z)；
S8.2)E-3DTV正则用以高光谱图像压缩重构的计算方法；
S8.2.1)E-3DTV正则用以高光谱图像压缩重构的优化模型如下
将Z分解成待修复高光谱干净数据项X和噪音项E,对X干净数据加E-3DTV正则，得到优化模型如下所示：
s.t.y＝Ψ(Z),Z＝X+E (17)
其中τ为调节参数，作用在于协调E-3DTV正则和误差项的作用折中。，将(17)式展开即得：
s.t.Y＝Ψ(Z)，Z＝X+E (18)
S8.2.2)给出公式(18)的增广拉格朗日函数
其中，Mi,Γ1,Γ2是拉格朗日乘子，μ是一个大于0的数；
S8.2.3)建立交替方向乘子法的迭代格式与终止条件：
μk+1＝ρμk (27)
迭代终止条件为：
S8.2.4)对问题(19)、(20)、(21)、(22)、(23)进行求解，给出迭代的具体算式；
S8.2.5)设置迭代的初始值为：E0＝0，由奇异值分解方法作用在Ψ*(y)上产生，其中Ψ*是Ψ的共轭算子；
S8.2.6)进行(19)-(27)的迭代运算，直到迭代满足终止条件。
传统的3DTV正则仅考虑了高光谱梯度图的结构稀疏性，并没有更深层次的对梯度图先验知识进行更深入的挖掘。而本发明所提的E-3DTV正则利用子空间学习的方法，对于梯度图的基进行稀疏约束，从而很好的刻画了高光谱梯度图结构的相关稀疏性，这一点是之前TV类方法从未考虑到的。首先，利用差分算子对高光谱图的三个维度进行差分操作，得到原高光谱图像的梯度图，空间上的梯度图可以反映出图片的边缘结构，光谱维度的梯度图可以反映出相邻谱段的相似性。通过算法的不断迭代，梯度图的基及其梯度图在基图下的表示矩阵和噪音矩阵均不断更新，同时基、系数矩阵及其噪音同时作用，不断反解出待修复的高光谱图。本发明通过对梯度图的认知，对3DTV正则进行本质的改进，从而在高光谱图像上大幅度的提升针对此类图像的修复效果。
附图说明
图1为本发明的流程图。
图2为本发明E-3DTV正则的阐释示意图。
图3为本发明在梯度图的基上进行修复能够带来鲁棒性效果的示意图。
图4为本发明基于E-3DTV正则在Indian Pines高光谱数据上band#220的去噪效果。其中(a)为原图，(b)为噪音，(c)为修复图。
图5为本发明基于E-3DTV正则在DC mall高光谱数据上band#120的去噪效果。其中(a)为原图，(b)为噪音，(c)为修复图。
图6为本发明基于E-3DTV正则在真实的Indian Pines高光谱数据上band#108的去噪效果。其中(a)为噪音图，(b)为修复图。
图7为本发明基于E-3DTV正则在真实的Urban高光谱数据上band#108的去噪效果。其中(a)为噪音图，(b)为修复图。
图8为本发明基于E-3DTV正则在真实的Terrian高光谱数据上band#207的去噪效果。其中(a)为噪音图，(b)为修复图。
图9为本发明基于E-3DTV正则在真实的Lowal高光谱数据上band#109的去噪效果。其中(a)为噪音图，(b)为修复图。
图10为本发明基于E-3DTV正则在采样率为5％的情况下DCmall高光谱数据上band#5的重构效果。其中(a)为原图，(b)为重构图。
图11为本发明基于E-3DTV正则在采样率为5％的情况下Urban高光谱数据上band#160的重构效果。其中(a)为原图，(b)为重构图。
图12为本发明基于E-3DTV正则在采样率为5％的情况下Lowal高光谱数据上band#160的重构效果。其中(a)为原图，(b)为重构图。
图13为本发明基于E-3DTV正则在采样率为5％的情况下Moffett高光谱数据上band#80的重构效果。其中(a)为原图，(b)为重构图。
具体实施方式
下面结合附图和实施例详细说明本发明的实施方式。
如图1所示，本发明一种基于E-3DTV正则的高光谱图像修复方法，包括如下步骤：
步骤S1：将含噪音的原三维高光谱数据沿着光谱维度展开成矩阵得到Y∈Rmn×T，其中m,n代表每个谱段下图像的长和宽，T代表高光谱的光谱数，初始化噪音项ε∈Rm×n×T的矩阵形式E∈Rmn×T和待修复的高光谱数据的矩阵表示X∈Rmn×T以及ADMM框架下的其他模型变量及参数，其他模型变量及参数包括ADMM中需要用到的上升乘子矩阵M和上升乘子常数的初始值μ和每次上升的倍率ρ，M＝0，μ的经验取值在[1e-3,1e-2]之间选择，ρ的经验取值在[1.05,1.8]之间选择；
步骤S2：对待修复的高光谱数据沿着水平、竖直、光谱三个维度做差分操作，得到的三个不同方向的梯度图，并沿着光谱维度展开成矩阵，记水平、竖直、光谱方向的梯度图分别为i＝1,2,3；
三个方向的差分图计算公式定义如下所示：
其中，i,j,k分别为高光谱图像对应的沿空间长度、宽度与谱段方向的坐标位置，表示原高光谱图像在位置(i,j,k)的灰度值，s＝1，2，3分别表示的是原高光谱的在位置的灰度值表示沿着空间长度、宽度与谱段对高光谱图进行差分运算之后得到的梯度图在位置(i，j，k)的取值。分别表示对应高光谱图像的(i,j,k)位置沿着空间长度、宽度与谱段方向分别平移一个位置的像素灰度值。为了保证差分图和待复原图尺寸一致便于计算，对差分图做zero-padding处理，将梯度图沿着光谱维度展开为矩阵，得到如下的公式：
其中unfold表示沿着光谱方向展开。
步骤S3：对三个方向的梯度图矩阵分别做低秩UV分解，其中基矩阵U反映了梯度图的共有的结构和特征信息，V是系数矩阵，它反映了梯度图矩阵的每一列再由基矩阵表达时的系数；对梯度图的基矩阵U加稀疏性约束，得到E-3DTV正则；
E-3DTV正则的形式如下：
其中，Ui为梯度矩阵的基，由于Ui中包含了梯度矩阵的共有结构或者特征，因此给Ui加上了一范数‖Ui‖1来约束基的稀疏性；Vi为梯度矩阵在基矩阵下的系数矩阵。为了保证上述优化问题解的唯一性，给系数矩阵加上ViVi T＝I的限制，这里I表示单位矩阵；r为低秩UV分解中的秩，r<min(mn,T)，蕴含了重建矩阵的低秩性信息。
由线性代数的知识可得，每一个谱段的梯度都可以通过基和系数来进行线性表示。相对传统的刻画梯度图稀疏性先验的3DTV而言，本发明的正则不是在上加一范数进行稀疏性约束，而是在基矩阵Ui上加一范数进行稀疏性约束。这样做的好处是可以充分的挖掘梯度图之间的相关性，通过稀疏基矩阵Ui和正交系数来表达梯度图的稀疏性和关联性，直观有效。
步骤S4：根据步骤S3，对待修复的X加E-3DTV正则，写出优化模型，利用ADMM框架进行迭代求解，到迭代稳定，获得已修复的图像X和噪音E。
优化模型如下所示：
s.t.Y＝X+E (4)
其中τ为调节参数，作用在于协调E-3DTV正则和误差项的作用折中，Y为输入的带噪音数据，E为噪音项，X为待复原的高光谱图数据，等式(4)展开即可得到下面的公式(5)：
利用ADMM框架进行迭代求解的步骤如下：
S6.1)给出公式(5)的增广拉格朗日函数：
其中，Mi,Γ是拉格朗日乘子，μ是上升乘子常数，<.,.>表示内积运算，表示矩阵的F范数；
S6.2)建立交替方向乘子法的迭代格式与终止条件：
μk+1＝ρμk (12)
迭代终止条件为：
S6.3)对问题(6)、(7)、(8)、(9)进行求解，给出迭代的具体算式；
S6.4)设置迭代的初始值为：E0＝0，由奇异值分解方法作用在Y上产生；
S6.5)进行公式(6)-(12)的迭代运算，直到迭代满足终止条件。
步骤S5：步骤S1到S4，给出了E-3DTV正则的一般高光谱图像修复流程，具体到去噪问题中直接套用。针对压缩重构问题，需要多加一步由少量采样数据y来计算重构数据的过程，具体地，高光谱图像压缩重构的具体流程如下：
S8.1)压缩重构问题简述
构造采样算子Ψ＝D.H.P，其中D是随机下采样算子，H是随机置换矩阵，P是WalshHadamard变换算子，固定采样算子Ψ之后，高光谱数据Z∈Rmn×T经Ψ的作用会得到采样过后的观测子y∈Rl,l＜＜mnT，压缩重构的目的就是通过观测y重构出Z。压缩重构蕴含内在等式约束为y＝Ψ(Z)；
S8.2)E-3DTV正则用以高光谱图像压缩重构的计算方法；
S8.2.1)E-3DTV正则用以高光谱图像压缩重构的优化模型如下
将Z分解成待复原高光谱数据项X和噪音项E,对X干净数据加E-3DTV正则，得到优化模型如下所示：
s.t.y＝Ψ(Z),Z＝X+E (17)
其中τ为调节参数，作用在于协调E-3DTV正则和误差项的作用折中。，将(17)式展开即得：
S8.2.2)给出公式(18)的增广拉格朗日函数
其中，Mi,Γ1,Γ2是拉格朗日乘子，μ是一个大于0的数；
S8.2.3)建立交替方向乘子法的迭代格式与终止条件：
μk+1＝ρμk (27)
迭代终止条件为：
S8.2.4)对问题(19)、(20)、(21)、(22)、(23)进行求解，给出迭代的具体算式；
1:更新噪音项E，从拉格朗日式子中，提取含E的元素，得到
2:更新基矩阵Ui,i＝1,2,3，从拉格朗日式子中，提取含Ui的项目，得到
3:更新表达系数Vi,i＝1,2,3，从拉格朗日式子中，提取含Vi的项目，得到
该问题可由矩阵的von Neumann’s trace inequality得到显式解:
4:更新X，从拉格朗日式子中，提取含X的项目，得到
对以上算式求导，将含X的式子移到等式一边，可得以下线性方程组。
与去噪模型相似，可算得其显式解形式，
5:更新Z，从拉格朗日式子中，提取含Z的项目，得到
对(33)求导，将含Z的式子移到等式一边，可得以下线性方程组
μΨ*ΨZ+μZ＝μ(X+E)+μΨ*y+Ψ*Γ1-Γ2 (34)
其中Ψ*是Ψ的伴随算子。该线性系统可以由现成的共轭梯度算法进行求解。
S8.2.5)设置迭代的初始值为：E0＝0，由奇异值分解方法作用在Ψ*(y)上产生，其中Ψ*是Ψ的共轭算子；
S8.2.6)进行(19)-(27)的迭代运算，直到迭代满足终止条件。
本发明的原理及效果可参考图2-图13。
图2为本发明E-3DTV正则的阐释示意图，图2中第a行的分别表示的是原高光谱图像数据及其沿着光谱维度展开后的矩阵图。第b行表示的是三个方向的差分操作示意图。第c行表示的是三个方向上得到的梯度图。第d行表示的是梯度图沿着光谱维度展开后得到的矩阵图。第e行表示的是梯度图的数值频数分布和矩阵的奇异值曲线，两条曲线表明梯度图具有很强的稀疏性和相关性，这就说明梯度图本身有很强的冗余性，这种冗余在这里就体现为相关性，而传统的3DTV只是刻画梯度图的稀疏性，并不能刻画这种相关性。第f行表示的对原图做了一个正交变化之后，可以得到一个更紧致的数据表达，这个表达实际上就是梯度矩阵的基。通过频数分布图可得知在梯度图的基上面也是具有稀疏性的。
在上面的描述，第f行得到的数据表达就是梯度图的基的证明如下：
对矩阵Gi做低秩矩阵分解，得到Gi＝UVT,其中VTV＝I。把U记为矩阵的基，V记为系数。等式两边同时乘以VT，可以得到GiV＝UVTV＝U。
图3为本发明在梯度图的基上面进行修复能够带来鲁棒性效果的示意图，可以阐明在梯度图的基上面进行修复能够带来的鲁棒性的效果。
图3中展示的第二列是干净图像的信息，第一列是加了少许噪音的图像。第a行是原始高光谱图像；第b行是水平方向的差分图信息，第c行是光谱维度上的差分图；通过第b,c两行可知在混杂噪音情形，图像的梯度被严重破坏；第d行是对水平方向高光谱图的梯度图做主成分提取得到的前两个基图；第e行是对光谱方向上的高光谱图的梯度图做主成分提取得到的前两个基图；通过第d,e两行可以知道在噪音情况，梯度图的基仍具有一定的鲁棒性。所以在构造E-3DTV的时候，选择在基图上进行操作是合理的。
图4为展示的E-3DTV方法在仿真实验Indian Pines数据上band#220的去噪的效果。第一列是原图，第二列是在原图中加上高斯核稀疏噪音得到的噪音图，第三列是E-3DTV方法去噪后的结果，该结果与原图越接近，说明去噪方法越好。为了方便展示本发明方法在高光谱图像上的修复效果，对图片的局部进行了放大。
图5展示的是E-3DTV方法在仿真实验DCmall数据上band#120的去噪的效果。第一列是原图，第二列是在原图中加上高斯核稀疏噪音得到的噪音图，第三列是E-3DTV方法去噪后的结果。同样为了方便展示本发明方法在高光谱图像上的修复效果，对图片的局部进行了放大。
图6展示的是E-3DTV方法在真实实验Indian Pines数据上band#108的去噪的效果。由于是真实数据，没有真实的无噪音的图片，所以只能通过修复后的视觉效果来进行观察。
图7展示的是E-3DTV方法在真实实验Urban数据上band#108的去噪的效果。同样为了方便展示本发明方法在高光谱图像上的修复效果，对图片的局部进行了放大。
图8和图9展示的是E-3DTV方法在真实实验Terrian数据上band#207和Lowal数据上band#109的去噪的效果。可以看到这两个数据集和前面的数据的结构以及真实的噪音是大不相同，而本发明方法依旧可以很好的去除噪音。同样为了方便展示本发明方法在高光谱图像上的修复效果，对图片的局部进行了放大。
图10到图13展示的是E-3DTV方法在高光谱压缩重构上的修复效果。同样为了展示方便，对修复图片的局部进行了放大。在展示的四幅图中，采样率均为5％，即本发明方法需要从原图的5％的采样数据修修复始图像。修复的图像和原始图像越相近，则修复的效果越好。其中图10展示的是DCmall数据上在band#5下的修复效果，图11展示的是Urban数据上在band#160下的修复效果，图12展示的是Lowal数据上在band#160下的修复效果，图12展示的是Moffett数据上在band#80下的修复效果。
Claims (9)
1.一种基于E-3DTV正则的高光谱图像修复方法，其特征在于，包括如下步骤：
步骤S1：将含噪音的原三维高光谱数据沿着光谱维度展开成矩阵得到Y∈Rmn×T，其中m,n代表每个谱段下图像的长和宽，T代表高光谱的光谱数，初始化噪音项ε∈Rm×n×T的矩阵形式E∈Rmn×T和待修复的高光谱数据的矩阵表示X∈Rmn×T以及ADMM框架下的其他模型变量及参数；
步骤S2：对待修复的高光谱数据沿着水平、竖直、光谱三个维度做差分操作，差分算子记作得到的三个不同方向的梯度图，并沿着光谱维度展开成矩阵，记水平、竖直、光谱方向的梯度图分别为
步骤S3：对三个方向的梯度图矩阵分别做低秩UV分解，即矩阵可以由两个低秩矩阵的乘积来近似表示，其中矩阵中基矩阵U反映了梯度图的共有结构和特征信息，V是系数矩阵，它反映了梯度图矩阵的每一列再由基矩阵表达时的系数；对梯度图的基矩阵U加稀疏性约束，得到E-3DTV正则；
步骤S4：根据步骤S3，对待修复的高光谱图像X加E-3DTV正则，写出优化模型，利用ADMM框架进行迭代求解，到迭代稳定，获得修复图像和对应噪音。
2.根据权利要求1所述基于E-3DTV正则的高光谱图像修复方法，其特征在于，所述步骤S1中，ADMM框架下的其他模型变量及参数包括ADMM中需要用到的上升乘子矩阵M和上升乘子常数的初始值μ和每次上升的倍率ρ，M＝0，μ的经验取值在[1e-3,1e-2]之间选择，ρ的经验取值在[1.05,1.8]之间选择。
3.根据权利要求1所述基于E-3DTV正则的高光谱图像修复方法，其特征在于：所述步骤S2中三个方向的差分图计算公式定义如下所示：
其中，i,j,k分别为高光谱图像对应的沿空间长度、宽度与谱段方向的坐标位置，表示原高光谱图像在位置(i,j,k)的灰度值， 分别表示沿着空间长度、宽度与谱段对高光谱图进行差分运算之后得到的梯度图在位置(i,j,k)的取值，分别表示对应高光谱图像的(i,j,k)位置沿着空间长度、宽度与谱段方向分别平移一个位置的像素灰度值，为了保证差分图和待修复图尺寸一致便于计算，对差分图做zero-padding处理，将梯度图沿着光谱维度展开为矩阵，得到如下的公式：
其中unfold表示沿着光谱方向展开。
4.根据权利要求1所述基于E-3DTV正则的高光谱图像修复方法，其特征在于，所述步骤S3中，E-3DTV正则的形式如下：
其中，Ui为梯度矩阵的基，由于Ui中包含了梯度矩阵的共有结构或者特征，因此给Ui加上了一范数‖Ui‖1来约束基的稀疏性；Ui为梯度矩阵在基矩阵下的系数矩阵，为了保证上述优化问题解的唯一性，给系数矩阵加上的限制，这里I表示单位矩阵；r为低秩UV分解中的秩，r<min(mn，T)，蕴含了重建矩阵的低秩性信息。
5.根据权利要求4所述基于E-3DTV正则的HSI图像修复方法，其特征在于，所述步骤S4中，建立的优化模型如下所示：
s.t.Y＝X+E (4)
其中τ为调节参数，用于协调E-3DTV正则和误差项的作用折中，‖E‖1表示的是对E加一范数，等式(4)展开即可得到下面的公式(5)：
6.根据权利要求5所述基于E-3DTV正则的高光谱图像修复方法，其特征在于，所述步骤S4中，利用ADMM框架进行迭代求解的步骤如下：
S6.1)给出公式(5)的增广拉格朗日函数：
其中，Mi，Γ是拉格朗日乘子，μ是上升乘子常数，<.,.>表示内积运算，表示矩阵的F范数；
S6.2)建立交替方向乘子法的迭代格式与终止条件：
μk+1＝ρμk (12)
迭代终止条件为：
S6.3)对问题(6)、(7)、(8)、(9)进行求解，给出迭代的具体算式；
S6.4)设置迭代的初始值为：E0＝0，由奇异值分解方法作用在Y上产生；
S6.5)进行公式(6)-(12)的迭代运算，直到迭代满足终止条件。
7.根据权利要求6所述基于E-3DTV正则的高光谱图像修复方法，其特征在于，所述步骤S6.3)中，每个问题的解如下所示：
1:更新待修复高光谱图X，从拉格朗日式子中，提取含X的元素，得到：
其中表示的是的共轭算子，fftn和ifftn分别为傅里叶变化和逆傅里叶变化，由于fftn是对三维张量操作的，所以方便运算起见增加了Fold算子，表示将矩阵沿着光谱拼回张量的操作；
2:利用软阈值收缩算子可得Ui显式解为：
其中软阈值收缩算子的定义如下所示：
其中阈值Δ为一个大于0的常数
3:由矩阵的奇异值分解和矩阵的迹不等式可得Vi显式解为：
4：由软阈值收缩算子可得到噪音项E：
8.根据权利要求1所述基于E-3DTV正则的高光谱图像修复方法，其特征在于，通过所述步骤S1-S4实现高光谱图像去噪；在所述步骤S1-S4的基础上，增加一步由少量采样数据y来计算重构数据的过程，可实现高光谱图像压缩重构。
9.根据权利要求8所述基于E-3DTV正则的高光谱图像修复方法，其特征在于，所述高光谱图像压缩重构的具体流程如下：
S9.1)压缩重构问题简述
构造采样算子Ψ＝D.H.P，其中D是随机下采样算子，H是随机置换矩阵，P是WalshHadamard变换算子，固定采样算子Ψ之后，高光谱原始数据Z∈Rmn×T经Ψ的作用会得到采样过后的观测子y∈Rl,l＜＜mnT，<<表示远远小于的符号压缩重构的目的就是通过观测y重构出Z，压缩重构蕴含内在等式约束为y＝Ψ(Z)；
S9.2)E-3DTV正则用以高光谱图像压缩重构的计算方法；
S9.2.1)E-3DTV正则用以高光谱图像压缩重构的优化模型如下
将Z分解成待修复高光谱干净数据项X和噪音项E,对X干净数据加E-3DTV正则，得到优化模型如下所示：
s.t.y＝Ψ(Z)，Z＝X+E (17)
其中τ为调节参数，用于协调E-3DTV正则和误差项的作用折中，‖E‖2表示的对噪音加二范数的约束，将(17)式展开即得：
S9.2.2)给出公式(18)的增广拉格朗日函数
其中，Mi,Γ1,Γ2是拉格朗日乘子，μ是一个大于0的数；
S9.2.3)建立交替方向乘子法的迭代格式与终止条件：
μk+1＝ρμk (27)
迭代终止条件为：
S9.2.4)对问题(19)、(20)、(21)、(22)、(23)进行求解，给出迭代的具体算式；
1:更新噪音项E，从拉格朗日式子中，提取含E的元素，得到
2:更新基矩阵Ui,i＝1,2,3，从拉格朗日式子中，提取含Ui的项目，得到
3:更新表达系数Vi,i＝1,2,3，从拉格朗日式子中，提取含Vi的项目，得到
该问题可由矩阵的von Neumann’s trace inequality得到显式解:
4:更新X，从拉格朗日式子中，提取含X的项目，得到
对以上算式求导，将含X的式子移到等式一边，可得以下线性方程组。
与去噪模型相似，可算得其显式解形式，
5:更新Z，从拉格朗日式子中，提取含Z的项目，得到
对(33)求导，将含Z的式子移到等式一边，可得以下线性方程组
μΨ*ΨZ+μZ＝μ(X+E)+μΨ*y+Ψ*Γ1-Γ2 (34)
其中Ψ*是Ψ的伴随算子。该线性系统可以由现成的共轭梯度算法进行求解。
S9.2.5)设置迭代的初始值为：E0＝0，由奇异值分解方法作用在Ψ*(y)上产生，其中Ψ*是Ψ的共轭算子；
S9.2.6)进行(19)-(27)的迭代运算，直到迭代满足终止条件。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201811046032.XA CN109191404B (zh) | 2018-09-07 | 2018-09-07 | 一种基于e-3dtv正则的高光谱图像修复方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201811046032.XA CN109191404B (zh) | 2018-09-07 | 2018-09-07 | 一种基于e-3dtv正则的高光谱图像修复方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN109191404A true CN109191404A (zh) | 2019-01-11 |
| CN109191404B CN109191404B (zh) | 2020-06-26 |
Family
ID=64915465
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201811046032.XA Active CN109191404B (zh) | 2018-09-07 | 2018-09-07 | 一种基于e-3dtv正则的高光谱图像修复方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN109191404B (zh) |
Cited By (10)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111127588A (zh) * | 2019-12-26 | 2020-05-08 | 中国人民解放军海军航空大学青岛校区 | 基于DirectX的大数据量参数曲线回放方法 |
| CN111223049A (zh) * | 2020-01-07 | 2020-06-02 | 武汉大学 | 一种基于结构-纹理分解的遥感图像变分融合方法 |
| CN111325697A (zh) * | 2020-03-04 | 2020-06-23 | 西安交通大学 | 一种基于张量本征变换的彩色图像修复方法 |
| CN111714124A (zh) * | 2020-06-18 | 2020-09-29 | 中国科学院深圳先进技术研究院 | 磁共振电影成像方法、装置、成像设备及存储介质 |
| CN111951182A (zh) * | 2020-07-14 | 2020-11-17 | 浙江工业大学 | 基于流形加权非局部曲率正则化的高光谱图像修复方法 |
| CN112634167A (zh) * | 2020-12-29 | 2021-04-09 | 南京理工大学 | 一种全变差协同范数约束迭代投影的高光谱图像滤波方法 |
| CN113160069A (zh) * | 2021-02-26 | 2021-07-23 | 桂林电子科技大学 | 一种基于图信号的高光谱图像去噪方法 |
| CN114972122A (zh) * | 2022-07-27 | 2022-08-30 | 中国科学院空天信息创新研究院 | 高光谱遥感图像坏像元复原方法、装置、电子设备及介质 |
| CN115034978A (zh) * | 2022-05-23 | 2022-09-09 | 南京邮电大学 | 基于子空间表示结合图嵌入约束的高光谱图像去噪方法 |
| CN116433534A (zh) * | 2023-06-09 | 2023-07-14 | 四川工程职业技术学院 | 一种高光谱图像修复方法、装置、存储介质及电子设备 |
Citations (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN108133465A (zh) * | 2017-12-29 | 2018-06-08 | 南京理工大学 | 基于空谱加权tv的非凸低秩松弛的高光谱图像恢复方法 |
-
2018
- 2018-09-07 CN CN201811046032.XA patent/CN109191404B/zh active Active
Patent Citations (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN108133465A (zh) * | 2017-12-29 | 2018-06-08 | 南京理工大学 | 基于空谱加权tv的非凸低秩松弛的高光谱图像恢复方法 |
Non-Patent Citations (2)
| Title |
|---|
| LE SUN 等: "Hyperspectral Restoration Employing Low Rank and 3D Total Variation Regularization", 《2016 INTERNATIONAL CONFERENCE ON PROGRESS IN INFORMATICS AND COMPUTING (PIC)》 * |
| 杨斌 等: "约束非负矩阵分解框架下高维自适应粒子群端元提取", 《遥感学报》 * |
Cited By (18)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111127588A (zh) * | 2019-12-26 | 2020-05-08 | 中国人民解放军海军航空大学青岛校区 | 基于DirectX的大数据量参数曲线回放方法 |
| CN111223049A (zh) * | 2020-01-07 | 2020-06-02 | 武汉大学 | 一种基于结构-纹理分解的遥感图像变分融合方法 |
| CN111223049B (zh) * | 2020-01-07 | 2021-10-22 | 武汉大学 | 一种基于结构-纹理分解的遥感图像变分融合方法 |
| CN111325697A (zh) * | 2020-03-04 | 2020-06-23 | 西安交通大学 | 一种基于张量本征变换的彩色图像修复方法 |
| CN111714124A (zh) * | 2020-06-18 | 2020-09-29 | 中国科学院深圳先进技术研究院 | 磁共振电影成像方法、装置、成像设备及存储介质 |
| CN111714124B (zh) * | 2020-06-18 | 2023-11-03 | 中国科学院深圳先进技术研究院 | 磁共振电影成像方法、装置、成像设备及存储介质 |
| CN111951182B (zh) * | 2020-07-14 | 2024-03-29 | 浙江工业大学 | 基于流形加权非局部曲率正则化的高光谱图像修复方法 |
| CN111951182A (zh) * | 2020-07-14 | 2020-11-17 | 浙江工业大学 | 基于流形加权非局部曲率正则化的高光谱图像修复方法 |
| CN112634167A (zh) * | 2020-12-29 | 2021-04-09 | 南京理工大学 | 一种全变差协同范数约束迭代投影的高光谱图像滤波方法 |
| CN112634167B (zh) * | 2020-12-29 | 2023-09-01 | 南京理工大学 | 一种全变差协同范数约束迭代投影的高光谱图像滤波方法 |
| CN113160069A (zh) * | 2021-02-26 | 2021-07-23 | 桂林电子科技大学 | 一种基于图信号的高光谱图像去噪方法 |
| CN113160069B (zh) * | 2021-02-26 | 2024-02-02 | 桂林电子科技大学 | 一种基于图信号的高光谱图像去噪方法 |
| CN115034978B (zh) * | 2022-05-23 | 2024-04-05 | 南京邮电大学 | 基于子空间表示结合图嵌入约束的高光谱图像去噪方法 |
| CN115034978A (zh) * | 2022-05-23 | 2022-09-09 | 南京邮电大学 | 基于子空间表示结合图嵌入约束的高光谱图像去噪方法 |
| CN114972122A (zh) * | 2022-07-27 | 2022-08-30 | 中国科学院空天信息创新研究院 | 高光谱遥感图像坏像元复原方法、装置、电子设备及介质 |
| CN114972122B (zh) * | 2022-07-27 | 2022-11-01 | 中国科学院空天信息创新研究院 | 高光谱遥感图像坏像元复原方法、装置、电子设备及介质 |
| CN116433534B (zh) * | 2023-06-09 | 2023-08-22 | 四川工程职业技术学院 | 一种高光谱图像修复方法、装置、存储介质及电子设备 |
| CN116433534A (zh) * | 2023-06-09 | 2023-07-14 | 四川工程职业技术学院 | 一种高光谱图像修复方法、装置、存储介质及电子设备 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN109191404B (zh) | 2020-06-26 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN109191404A (zh) | 一种基于e-3dtv正则的高光谱图像修复方法 | |
| CN111260576B (zh) | 一种基于去噪三维卷积自编码网络的高光谱解混算法 | |
| CN109035142B (zh) | 一种对抗网络结合航拍图像先验的卫星图像超分辨方法 | |
| CN111369487B (zh) | 一种高光谱和多光谱图像融合方法、系统及介质 | |
| US20130336540A1 (en) | Decomposition apparatus and method for refining composition of mixed pixels in remote sensing images | |
| CN106709881A (zh) | 一种基于非凸低秩矩阵分解的高光谱图像去噪方法 | |
| CN105931264B (zh) | 一种海面红外小目标检测方法 | |
| CN110400276B (zh) | 高光谱图像去噪方法、装置 | |
| CN106447630B (zh) | 基于概率矩阵分解的高光谱图像锐化方法 | |
| CN103279959B (zh) | 一种二维分析稀疏模型、其字典训练方法和图像去噪方法 | |
| CN112712034B (zh) | 一种高光谱图像的解混方法及系统 | |
| CN110139046B (zh) | 一种基于张量的视频帧合成方法 | |
| CN103020939A (zh) | 利用多时相数据去除光学遥感影像大面积厚云的方法 | |
| Cao et al. | CS-MRI reconstruction based on analysis dictionary learning and manifold structure regularization | |
| CN103413351B (zh) | 基于压缩感知理论的三维人脸快速重建方法 | |
| CN106296583B (zh) | 基于图像块组稀疏编码与成对映射的含噪高光谱图像超分辨率重构方法 | |
| CN111292266A (zh) | 基于双低秩矩阵分解的gf-5遥感影像混合噪声去除方法 | |
| Hesabi et al. | Structure and texture image inpainting | |
| CN116626765A (zh) | 一种基于二维k-svd和卷积稀疏编码的多道地震反褶积方法 | |
| CN113094645B (zh) | 一种形态成分约束优化的高光谱数据解混方法 | |
| Karoui et al. | A Novel Linear Mixing Model Addressing Spectral-Spatial Intra-Class Variability with an Associated Penalized NMF-Based Hyperspectral Unmixing Algorithm | |
| CN108460777A (zh) | 一种面向植物高光谱的提取分块压缩重构方法 | |
| CN112634167A (zh) | 一种全变差协同范数约束迭代投影的高光谱图像滤波方法 | |
| CN114359064B (zh) | 基于双重梯度约束的高光谱图像恢复方法 | |
| CN116757925B (zh) | 一种生成高时空谱分辨率卫星遥感影像的方法和装置 |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |