[发明专利]一种变分PDE图像保边保结构复原方法

权利要求书

1.一种变分PDE图像保边保结构复原方法，包括以下步骤：

(1)输入一幅大小为N×N的退化图像f(x,y)，令f(i,j)表示图像在第i行第j列的像素值，利用Matlab中的diff函数定义函数g(x,y)的梯度算子和海森算子

▿2g=(▿xxg,▿xyg;▿yxg,▿yyg)T:]]>

▿xg=[diff(g,1,2),g(:,1)-g(:,N)];]]>

▿yg=[diff(g,1,1);g(1,:)-g(m,:)];]]>

▿xxg=[▿xg(:,1)-▿xg(:,N);diff(▿xg,1,2)];]]>

▿yyg=[▿yg(1,:)-▿yg(N,:);diff(▿yg,1,1)];]]>

▿xyg=[diff(▿xg,1,1);▿xg(1,:)-▿xg(N,:)].]]>

选取初始值u⁰＝f,w⁰＝0,v⁰＝0,z⁰＝0,α⁰＝0,β⁰＝0,λ⁰＝0；

(2)求解离散化模型

其中

|▿u|l1=Σi=1NΣj=1N|▿ui,j|=Σi=1NΣj=1N(▿xui,j)2+(▿yui,j)2,]]>

|▿2u|l1=Σi=1NΣj=1N|▿2ui,j|=Σi=1NΣj=1N(▿xxui,j)2+(▿xyui,j)2+(▿yxui,j)2+(▿yyui,j)2;]]>

即引入约束条件

w=u,v=▿u,z=▿2u;]]>

将其转化为下面极大极小的鞍点问题

其中

这里α,β,λ分别为对应的拉格朗日乘子，c₁,c₂,c₃分别为对应的正的罚参数，为示性函数，即则否则，

(3)极大极小的鞍点问题中耦合了七个变量，在算法中应用Gauss-Seidl迭代格式，

求解子问题u^k+1，即考虑

由于w^k,v^k,z^k,α^k,β^k,λ^k可看做常变量，因此舍弃与u无关的变量便有

求出其对应的显式解，即对应的Euler-Lagrangian方程为

其中div²为的伴随矩阵，由(1)中选取的为循环边界差分格式，同时利用散度定理中的和因此可以利用快速傅里叶变换求解u^k+1，即为：其中表示快速傅里叶变换，表示快速傅里叶变换的逆变换；上述算子具体实现在程序中定义为：

conjoI＝conj(psf2otf([1,0],[N,N]))；

conjoDx＝conj(psf2otf([1,-1],[N,N]))；

conjoDy＝conj(psf2otf([1；-1],[N,N]))；

conjoDxx＝conj(psf2otf([1,-2,1],[N,N]))；

conjoDyy＝conj(psf2otf([1；-2；1],[N,N]))；

conjoDxy＝conj(psf2otf([1,-1；-1,1],[N,N]))；

otfH＝psf2otf(K,[N,N])；

求解子问题w^k+1，即考虑

舍去与w无关项，则对应的子问题为

wk+1=argminw<α,w>l2c12||w-uk+1||l22+δC(w).]]>

该问题用投影梯度法快速求解，因此有

求解子问题v^k+1，即考虑

舍去与v无关项，则对应的子问题为

vk+1=argminvc22||v-(▿uk+1-βkc2)||l22+b|v|l1.]]>

该问题对应经典的压缩阈值问题，因此其显式解为

vk+1=▿uk+1-βkc2|▿uk+1-βkc2|l1max{|▿uk+1βkc2|l1-bc2,0}]]>

求解子问题z^k+1，由于子问题z与v具有类似的形式，因此同样可写成其对应的显式解为

zk+1=▿2uk+1-λkc3|▿2uk+1-λkc3|l1max{|▿2uk+1-λkc3|l1-1-bc3,0}]]>

由于拉格朗日乘子对应的问题是平凡的，因此其更新变量即为

α^k+1＝α^k+c₁(w^k+1-u^k+1)；

βk+1=βk+c2(vk+1-▿uk+1);]]>

λk+1=λk+c3(zk+1-▿2uk+1);]]>

终止标准，由交替方向乘子法的收敛性定理可知，当收敛到零时，算法的解收敛到初始问题的解，因此在算法中选取终止标准定为当算法终止时，输出复原图像u:＝u^k+1。