Replicated边界条件下图像模糊矩阵与图像矢量乘积的一种替代计算方法技术

本发明专利技术公开在Replicated边界条件下，图像模糊矩阵与图像矢量乘积的一种替代计算方法，包括：1)在Replicated边界条件下，将图像模糊矩阵与图像矢量乘积(乘积一)，以及图像模糊矩阵的转置与图像矢量乘积(乘积二)，分别化为多个带可利用结构的分块矩阵与图像矢量的乘积之和；2)对各个分块矩阵，构造对应的点扩展函数；3)用分块矩阵的点扩展函数与图像或图像对应边界间的卷积，替代地计算两种乘积的矩阵形式中的各个分解部分；以及4)计算两种乘积。所公开的计算方法，可集成到Replicated边界条件下的大型图像滤波与图像恢复中，用于解决其中图像模糊矩阵及其转置与图像矢量乘积难于计算的问题。

【技术实现步骤摘要】
Replicated边界条件下图像模糊矩阵与图像矢量乘积的一种替代计算方法
本专利技术涉及图像处理，特别地，本专利技术涉及到Replicated边界条件下的线性图像恢复中，对图像模糊矩阵及其转置与图像矢量乘积的一种替代计算。
技术介绍
图像恢复是对线性位移不变引起的下述退化问题进行求解f＝Hx+η(1)其中：图像模糊矩阵H∈Rmn×mn由点扩展函数(PointSpreadFunction,PSF)h∈Rp×q的元素构成，其结构唯一地由所采用的边界条件决定；f∈Rmn×1由退化图像F∈Rm×n按行优先排列得到的矢量；η为加性噪声的重排矢量。在对(1)的求解中，主要的计算任务,是计算H与f乘积(乘积一)，即g1＝Hf，以及H的伴随矩阵H*(在实数域中，即为H的转置矩阵HT)与f的乘积(乘积二)，即g2＝HTf。一般地，H和HT皆为严重稀疏的超大型矩阵，使得乘积一和乘积二难于直接计算。对此，针对一些边界条件下的乘积一和乘积二，目前已有两种替代计算方式：一是基于Kronecker乘积的替代计算方式，一种是基于对图像模糊矩阵进行对角化的替代计算方式。在确定性(Determinate)边界条件下，乘积一和乘积二能否被替代计算，以及能的话，被上述哪一种方式替代计算更为有效，这两个问题都取决于采用的是何种边界条件，以及点扩展函数是否对称。目前，确定性的图像边界条件有Zero边界条件、Periodic边界条件、Neumann边界条件、Anti-Reflective边界条件、Mean边界条件和Replicated边界条件，在这些条件下，图像模糊矩阵具有与退化图像不相关的确定性结构，且结构只唯一地由边界条件决定，这使得我们可以去分析这些结构，并将其应用到乘积一和乘积的计算,和对(1)式的预置器的设计中。Replicated边界条件是一种最新的边界条件，由我国学者Xuzhou等于2014年在文献”Aboundaryconditionbaseddeconvolutionframeworkforimagedeblurring”中提出，Xuzhou等认为该边界条件考虑了图像像素符合“重尾”分布的情况(即图像中大多像素和它的邻域像素值相近)，因此是一种较为符合图像实际的边界条件，可应用在对具有相应特征的图像的恢复和处理中。Xuzhou等在上述同一文献中还提出了基于Replicated边界的图像去模糊框架：1)对H进行乘法分解，即将H转化为非方的有偏模糊矩阵K∈R(m+2p1)(n+2q1)×mn与非方的边界条件矩阵P∈Rmn×(m+2p1)(n+2q1)的乘积，即H＝KP，从而将图像退化的问题f＝Hx+ζ转化为f＝KPx+ζ；2)采用共轭梯度算法求解线性方程P*K*y＝P*K*KPx，来恢复未知图像x，并在对应的迭代求解中，利用两次快速傅立叶变换(FastFourierTransformation，FFT)和一次快速傅立叶逆变换来计算乘积一和乘积二的矩阵形式；3)确定了Replicated边界条件下P和P*的结构，使得Px和P*x可被计算，从而提出了Replicated边界条件下乘积一和乘积二的一种替代计算方法，也是该边界条件下唯一的一种替代计算方法。XuZhou等对乘积一和乘积二的替代计算方法算法复杂度为O(m×n×Nnz)，其中Nnz为点扩展函数中非零元素的个数。在其他边界条件下，该框架是否可行，则取决于边界矩阵P和P*的结构是否已确定。
技术实现思路
本专利技术所要解决的技术问题是：设计Replicated边界条件下乘积一和乘积二的一种新的间接计算方法，该方法不对模糊矩阵进行乘法分解，而是进行加法分解，使对乘积一和乘积二的计算精度更高，并且在点扩展函数为常规大小时，时间耗费与XuZhou的方法相当。本专利技术针对上述问题，给出一种Replicated边界条件下乘积一和乘积二的间接计算方法，思路是：在Replicated边界条件下，根据图像模糊矩阵的结构，先将乘积一和乘积二，分别分解为多个矩阵与矢量乘积的和，保证所产生的分解矩阵带可利用的分块结构；利用各分块矩阵的结构，构造其卷积核；用卷积替代地计算各分解矩阵与图像矢量的乘积的矩阵形式；综合各计算结果，得到乘积一和乘积二。为支持上述思路，本专利技术中：实施例1给出Replicated边界条件下的乘积一和乘积二的替代计算方法的具体步骤；实施例2将Replicated边界条件下的乘积一和乘积二，分别转化为多个模糊矩阵与图像矢量的乘积之和；实施例3针对乘积一与乘积二的分解式中的各分解矩阵，构造对应的点扩展；实施例4基于构造的点扩展函数,采用卷积，计算乘积一和乘积二各个分解矩阵与图像矢量乘积的矩阵形式；实施例5通过实验结果，对本专利技术给出的卷积替代方法的有效性进行了说明。本专利技术所给出的计算方法，提供了Replicated边界条件下，对应大型乘积一和乘积二的替代计算，可集成到Replicated边界条件下的图像滤波和图像恢复中，具有应用价值。本专利技术所给出的乘积一和乘积二的计算精度较高。本算法的时间复杂度，主要来自于对零边界条件下的乘积一和乘积二的卷积运算复杂度，即O(m×n×p×q)，因此当点扩展函数为普通情形，即普通大小且元素值随机时，本方法与XuZhou的方法时间耗费相当，而当点扩展的大小较大而非零元素个数较少时，本专利方法会比XuZhou的方法慢。附图说明图1示出的本专利技术的计算流程；图2示出的是实验中所用的灰度图像Saturno图像，大小为512×512；图3示出的是实验中所用的灰度图像Cameraman图像,大小为256×256；图4示出是实验中图像Saturno所用的点扩展函数Saturnopsf,图5、图6和图7示出的是实验中图像Cameraman所用的三个随机生成的点扩展函数PSFn1,PSFn2和PSFn3,大小分别为53×53,7×7，17×13和31×37。具体实施方式在以下描述中，出于对专利说明的目的，结合附图，对
技术实现思路
中所述实施例的具体内容进行阐述。为了表达清晰和方便于编程实现，一些计算公式中应用了Matlab的函数形式。第1实施例如图1所示，本专利技术给出Replicated边界条件下乘积一g1＝Hf和乘积二g2＝HTf的一种替代计算方法，包括以下步骤：(1)输入F和h，对F按行优先重排成图像矢量f；(2)对g1和g2进行分解；(3)对g1和g2分解公式中各分解部分涉及的模糊矩阵，构造其点扩展函数，据此，采用点扩展函数与图像或对应图像边界间卷积，计算各分解部分的矩阵形式；(4)计算g1和g2的对应图像；(5)对g1和g2的对应图像，按行优先重排成矢量，得到g1和g2。第2实施例本专利技术对g1的分解为g1＝Hf＝(HTT+HTL+HLT+HLL)f,对g2的分解为g2＝HTf＝(HTTT+HTTL+HTLT+HTLL)f，两个分解公式中：各矩阵中，上标T表示对矩阵的转置，第1个下标代表当前矩阵的块间结构，第2个下标代表当前矩阵的块内结构，T和L分别具体代表Toeplitz矩阵和L-type矩阵，这里不给出T与L矩阵的具体结构，这不会影响后续对计算方法的描述。第3实施例本专利技术给出乘积一分解公式g1＝Hf＝(HTT+HTL+HLT+HLL)f和乘积二分解公式g2＝HTf＝(HTTT+H本文档来自技高网...

【技术保护点】
1.Replicated边界条件下图像模糊矩阵与图像矢量乘积的一种替代计算方法，所述方法在退化图像为F∈Rm×n和点扩展函数为h∈Rp×q的图像恢复中，当图像边界条件为Replicated边界条件时，对由h的元素构成的图像模糊矩阵H∈Rmn×mn与F的重排图像矢量f∈Rmn×1间的乘积，即乘积一g1＝Hf，以及H的转置HT与f间的乘积，即乘积二g2＝HTf，采用如下步骤替代地进行计算：(1)输入F和h，对F按行优先重排成图像矢量f；(2)对g1和g2计算公式中的模糊矩阵及其转置进行加法分解；(3)对g1和g2分解公式中各分解部分所涉及的矩阵，构造其点扩展函数，由各点扩展函数与图像或对应图像边界间的卷积，计算各分解部分的矩阵形式；(4)计算g1和g2的对应图像；(5)对g1和g2的对应图像，按行优先重排成矢量，得到g1和g2。

【技术特征摘要】
1.Replicated边界条件下图像模糊矩阵与图像矢量乘积的一种替代计算方法，所述方法在退化图像为F∈Rm×n和点扩展函数为h∈Rp×q的图像恢复中，当图像边界条件为Replicated边界条件时，对由h的元素构成的图像模糊矩阵H∈Rmn×mn与F的重排图像矢量f∈Rmn×1间的乘积，即乘积一g1＝Hf，以及H的转置HT与f间的乘积，即乘积二g2＝HTf，采用如下步骤替代地进行计算：(1)输入F和h，对F按行优先重排成图像矢量f；(2)对g1和g2计算公式中的模糊矩阵及其转置进行加法分解；(3)对g1和g2分解公式中各分解部分所涉及的矩阵，构造其点扩展函数，由各点扩展函数与图像或对应图像边界间的卷积，计算各分解部分的矩阵形式；(4)计算g1和g2的对应图像；(5)对g1和g2的对应图像，按行优先重排成矢量，得到g1和g2。2.根据权利要求1所述的图像处理方法，其特征在于对所述的对g1和g2进行分解，具体分解方法为：对g1的分解为g1＝Hf＝(HTT+HTL+HLT+HLL)f,对g2的分解为g2＝HTf＝(HTTT+HTTL+HTLT+HTLL)f，其中：各矩阵中，第1个下标代表当前矩阵的块间结构，第2个下标代表当前矩阵的块内结构，下标T和L分别具体代表Toeplitz矩阵和L-type矩阵，上标T表示对矩阵的转置。3.根据权利要求2所述的图像处理方法，其特征在于对g1和g2分解公式中,各分解部分所对应的矩阵，其点扩展函数的构造方法为：用fliplr(X)表示对X左右置反，flipud(X)表示对矩阵上下置反，V(:,:)中的:表示取到矩阵V的所有行或列，记和则：(1)HTT的点扩展函数为h1＝h；(2)构造HTL的点扩展函数，共2个核，方法为和其中i＝1,...,q1j＝1,...,p；(3)构造HLT的点扩展函数，共2个核，方法为和其中1≤i≤p1；(4)构造HLL的点扩展函数，共4个核，方法为和其中i＝1,...,p1；(5)HTTT的点扩展函数为h2＝fliplr(flipud(h))；(6)构造HTTL的点扩展函数，共2个核，方法为h2_TL_1＝fliplr((h1_TL_1)T)和h2_TL_2＝fliplr((h1_TL_2)T)；(7)构造HTLT的点扩展函数，共2个核，方法为h2_LT_1＝fliplr(flipud(h1_LT_1))和h2_LT_2＝fliplr(flipud(h1_LT_2))；(8)构造HTLL的点扩展函数，共4个核,方法为h2_LL_1＝fliplr(flipud((h1_LL_1)T))，h2_LL_2＝fliplr(flipud((h1_LL_2)T))，h2_LL_3＝fliplr(flipud((h1_LL_3)T))和h2_LL_4＝fliplr(flipud((h1_LL_4)T))。4.根据权利要求3所述的图像处理方法，其特征在于对g1和g2分解公式中的各分解部分，其矩阵形式的具体计算方法为：用conv2(.,.)表示二维卷积，用i:j表示从i到j的所有整数，用V(end,end)中的end标引对应矩阵V的最后行与最后列，用G1_TT、G1_TL、G1_LT、G1_LL、G2_TT、G2_TL、G2_LT和G2_LL，分别表示HTTf、HTLf、HLTf、HLLf、HTTTf、HTTLf、HTLTf和HTLLf等矩阵矢量乘积的矩阵形式，则：(1)G1_TT的计算方法为G1_TT＝conv2(F,h1,'same')；(2)G1_TL的计算方法为：计算r1_TL...

【专利技术属性】
技术研发人员：王令群，袁小华，王振华，常英立，李丛，张天蛟，查铖，邢鹏，谢沚娟，
申请(专利权)人：上海海洋大学，
类型：发明
国别省市：上海,31