CN105957024A - Blind deblurring method based on image block prior and sparse norm - Google Patents

Blind deblurring method based on image block prior and sparse norm Download PDF

Info

Publication number
CN105957024A
CN105957024A CN201610248012.5A CN201610248012A CN105957024A CN 105957024 A CN105957024 A CN 105957024A CN 201610248012 A CN201610248012 A CN 201610248012A CN 105957024 A CN105957024 A CN 105957024A
Authority
CN
China
Prior art keywords
image
candidate image
blur kernel
pyramid
block
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
Application number
CN201610248012.5A
Other languages
Chinese (zh)
Other versions
CN105957024B (en
Inventor
李阳阳
梁晓旭
王哲
焦李成
刘芳
尚荣华
马晶晶
刘若辰
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Xidian University
Original Assignee
Xidian University
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Xidian University filed Critical Xidian University
Priority to CN201610248012.5A priority Critical patent/CN105957024B/en
Publication of CN105957024A publication Critical patent/CN105957024A/en
Application granted granted Critical
Publication of CN105957024B publication Critical patent/CN105957024B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00Image enhancement or restoration
    • G06T5/73Deblurring; Sharpening
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/20Special algorithmic details
    • G06T2207/20021Dividing image into blocks, subimages or windows

Landscapes

  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • Image Analysis (AREA)

Abstract

本发明提出一种基于图像块先验与稀疏范数的盲去模糊方法,主要是解决现有技术对图像去模糊质量差的问题,其方案是:输入模糊图像;初始化模糊核、二进制掩模,候选图像;调用金字塔模型,将候选图像根据金字塔层数下采样,对候选图像与模糊核上采样;更新二进制掩模、更新图像块方差、更新图像样例块;固定参数更新模糊核、候选图像,直到金字塔最后一层;设置迭代次数,固定模糊核以及候选图像的范数保持不变,对模糊核添加的l1范数正则,得到新候选图像;固定候选图像不变,对候选图像添加的l1/l2范数正则项,得到新模糊核;直到迭代到最高次数。本发明提升了盲去模糊的效果和鲁棒性,可用于医疗器械、计算机视觉及图像视频处理。

The present invention proposes a blind deblurring method based on image block prior and sparse norm, which mainly solves the problem of poor image deblurring quality in the prior art. The solution is: input blurred image; initialize blur kernel and binary mask , the candidate image; call the pyramid model, downsample the candidate image according to the number of pyramid layers, upsample the candidate image and the blur kernel; update the binary mask, update the variance of the image block, update the image sample block; update the blur kernel, candidate Image, until the last layer of the pyramid; set the number of iterations, fix the blur kernel and the norm of the candidate image remain unchanged, and add l 1 norm regularization to the blur kernel to get a new candidate image; the fixed candidate image remains unchanged, and the candidate image Added l 1 /l 2 norm regularization term to get a new fuzzy kernel; until the highest number of iterations is reached. The invention improves the effect and robustness of blind deblurring, and can be used in medical equipment, computer vision and image and video processing.

Description

基于图像块先验与稀疏范数的盲去模糊方法Blind Deblurring Method Based on Image Block Prior and Sparse Norm

技术领域technical field

本发明属于图像处理技术领域,特别涉及一种图像盲去模糊方法,可用于航空航天、医疗器械、计算机视觉以及图像视频处理。The invention belongs to the technical field of image processing, in particular to an image blind deblurring method, which can be used in aerospace, medical equipment, computer vision, and image and video processing.

背景技术Background technique

随着个人智能设备的广泛推广使用,人手一台照相机已然成为事实。在这个智能化的时代,人们无时无刻不在与图片打交道,小到生活中自己拍摄的生活图片,大到神舟七号伴飞小卫星对神州七号的长达20分钟的拍照。这其中都会涉及到图像去模糊的相关知识,因此研究图像去模糊也就变得有相当重要的有意义。With the widespread use of personal smart devices, one camera per person has become a reality. In this era of intelligence, people deal with pictures all the time, ranging from the small pictures of their own life in their daily lives to the 20-minute photos taken by the small satellite that accompanied Shenzhou 7 to Shenzhou 7. This will involve the relevant knowledge of image deblurring, so it is very important and meaningful to study image deblurring.

图像在形成,记录,传输的过程中,由于受光学成像系统的相差,成像衍射,成像非线性,系统噪声多种因素的影响,图像的质量都会有所下降,图像的如上所述的一系列的过程即为图像的退化。而图像恢复,亦称为图像复原,就是尽最大可能减少或消除图像的质量的下降,恢复被退化图像的本来面目。In the process of image formation, recording and transmission, due to the influence of various factors such as phase difference of the optical imaging system, imaging diffraction, imaging nonlinearity and system noise, the quality of the image will decrease. The process is the degradation of the image. Image restoration, also known as image restoration, is to reduce or eliminate the degradation of image quality as much as possible, and restore the original appearance of degraded images.

初期的图像去模糊问题可以追溯至上世纪探索外太空兴起时期,由于受到极端恶劣天气或者大气湍流等因素影响,会造成图像质量的下降。这种图像质量的下降对于科学研究有非常大的影响。因此,图像去模糊的研究是非常有必要且又富有挑战性的研究内容。再比如在视频监控方面,监控视频的获得是不可逆的,也就是调阅之前的监控视频发现有模糊的场景存在时,监控视频是不能重新拍摄的,这时就需要进行图像去模糊才能获取更多的信息。此外,图像去模糊还在其他许多方面都有很重要的应用,比如材料科学图像处理,公安、历史、人文照片图像复原,扫描文档处理,星载,机载等航空侦察系统等方面。从单幅模糊图像中恢复出一幅清晰的图像早已成为一个日益基础且重要的研究问题。The initial image deblurring problem can be traced back to the rise of the exploration of outer space in the last century. Due to factors such as extreme weather or atmospheric turbulence, the image quality will be reduced. This degradation of image quality has a very large impact on scientific research. Therefore, the research on image deblurring is very necessary and challenging. Another example is in video surveillance, the acquisition of surveillance video is irreversible, that is, when a blurred scene is found in the previous surveillance video, the surveillance video cannot be re-shot. At this time, image deblurring is required to obtain more much information. In addition, image deblurring is also very important in many other aspects, such as material science image processing, public security, history, cultural photo image restoration, scanned document processing, spaceborne, airborne and other aerial reconnaissance systems, etc. Recovering a clear image from a single blurred image has become an increasingly fundamental and important research problem.

图像去模糊是与图像退化完全相反的过程,根据退化的成因是否已知,可以将图像去模糊任务分为图像非盲去模糊和图像盲去模糊两种。图像的非盲去模糊即是知道图像的退化原因,从而从一幅模糊图像中恢复出一幅清晰图像,通常需要注意的问题是减少可能存在的振铃效应和不可抑制的噪声。而图像盲去模糊则是在并不知道退化原因的情况下,从一幅模糊图像中恢复出一幅清晰的图像。图像的盲去模糊面临的问题更多,因为图像的退化原因以及清晰图像都是未知的。因此,对于图像去模糊问题的研究是非常有意义而且是非常有必要的。Image deblurring is a completely opposite process to image degradation. According to whether the cause of degradation is known, image deblurring tasks can be divided into two types: non-blind image deblurring and image blind deblurring. The non-blind deblurring of the image is to know the cause of image degradation, so as to restore a clear image from a blurred image. Usually, the problem that needs attention is to reduce the possible ringing effect and unsuppressible noise. Blind image deblurring is to restore a clear image from a blurred image without knowing the cause of the degradation. Blind deblurring of images is more problematic because both the cause of image degradation and the clear image are unknown. Therefore, the research on image deblurring is very meaningful and necessary.

图像模糊在日常生活中十分常见,对于模糊核已知的去模糊过程,称之为非盲去模糊;对于模糊核未知的去模糊过程,称之为盲去模糊。对整幅图像的每个像素都通过相同模糊方式进行模糊,称之为均匀模糊;反之在图像的不同区域,模糊的方式各不相同,称之为非均匀模糊。非均匀模糊更贴近于现实生活中的模糊现象。Image blurring is very common in daily life. For the deblurring process with known blur kernel, it is called non-blind deblurring; for the deblurring process with unknown blur kernel, it is called blind deblurring. Each pixel of the entire image is blurred by the same blurring method, which is called uniform blurring; on the contrary, in different areas of the image, the blurring methods are different, which is called non-uniform blurring. Non-uniform blur is closer to the blur phenomenon in real life.

对于均匀盲去模糊的求解过程,目标就是从模糊图像中恢复清晰图像和模糊核。显然这是一个病态问题,即有多组候选解。为了克服盲去模糊过程中的病态问题,早期的去模糊理论提出了参数化模糊核的方法,认为模糊核的形状为线性核,只有长度和角度两个变量,虽然该方法解决了一部分去模糊的问题,但是现实生活中的模糊成因相对比较复杂,所以存在相当大的应用局限性。之后有些前辈的工作中对模糊图像添加假设或者先验知识。例如,通常假设它是稀疏的并且是连续的,通常假设图像的梯度信息服从重尾分布。在此之后,一些新的理论方法取得了更好的效果,可以估计出更可靠的模糊核、图像质量更高的候选图像以及更好的鲁棒性。而后在商业软件上也出现了很好的应用,例如Photoshop系列软件提供的防抖功能。For the solution process of uniform blind deblurring, the goal is to recover the sharp image and the blur kernel from the blurred image. Obviously this is an ill-conditioned problem, that is, there are multiple sets of candidate solutions. In order to overcome the ill-conditioned problem in the blind deblurring process, the early deblurring theory proposed a method of parameterizing the blurring kernel. It believed that the shape of the blurring kernel is a linear kernel with only two variables, length and angle. Although this method solves part of the deblurring However, the causes of fuzzy in real life are relatively complex, so there are considerable application limitations. Later, some predecessors added assumptions or prior knowledge to blurred images in their work. For example, it is usually assumed that it is sparse and continuous, and it is usually assumed that the gradient information of the image obeys a heavy-tailed distribution. After this, some new theoretical methods achieved better results, estimating more reliable blur kernels, candidates with higher image quality, and better robustness. Then there are good applications in commercial software, such as the anti-shake function provided by Photoshop series software.

目前,盲去模糊方法进入了一个发展的黄金期,主要有三种方法:At present, the blind deblurring method has entered a golden period of development, and there are three main methods:

第一种为基于最大后验概率的方法,该类方法寻求最有可能的解决方案,最大限度地提高模糊核k和候选图像x的联合后验概率分布。该类方法简单易懂,缺点是有时可能会收敛到我们并不希望的解;The first is the method based on the maximum posterior probability, which seeks the most likely solution to maximize the joint posterior probability distribution of the blur kernel k and the candidate image x. This type of method is simple and easy to understand, but the disadvantage is that it may sometimes converge to a solution we do not want;

第二种为基于变分贝叶斯的方法,该类方法是在基于最大后验概率的基础上发展出来的,由于其考虑所有可能的解,因而比基于最大后验概率的算法鲁棒性更好,但缺点是速度较慢;The second is the method based on variational Bayesian, which is developed on the basis of the maximum posterior probability. Because it considers all possible solutions, it is more robust than the algorithm based on the maximum posterior probability. better, but at the disadvantage of being slower;

第三种是基于边缘预测的方法,该类方法认为模糊核k可以从一小部分图像边缘中估计出来,并使用启发式图像滤波器来恢复锐利的边缘,该类方法在模糊核k估计阶段速度非常快,并且在实验中证明是有效的,但是因为加入了启发式滤波的步骤,所以理论分析特别困难。The third method is based on edge prediction. This method considers that the blur kernel k can be estimated from a small part of the image edge, and uses a heuristic image filter to restore the sharp edge. This type of method is in the blur kernel k estimation stage The speed is very fast, and it is proved to be effective in the experiment, but because the step of heuristic filtering is added, the theoretical analysis is particularly difficult.

发明内容Contents of the invention

本发明针对上述方法中的不足,提出一种基于图像块先验与稀疏范数的盲去模糊方法以提升盲去模糊的适应性、可靠性和鲁棒性。Aiming at the deficiencies in the above methods, the present invention proposes a blind deblurring method based on image block prior and sparse norm to improve the adaptability, reliability and robustness of blind deblurring.

本发明的技术关键是:在金字塔模型的最后一层中,添加一项基于当下效果最好的稀疏范数正则项,使得迭代过程中,估计模糊核的阶段指向正确的方向,从而得到更加逼近真实场景的模糊核,其实现步骤包括如下:The technical key of the present invention is: in the last layer of the pyramid model, add a sparse norm regularization item based on the best current effect, so that in the iterative process, the stage of estimating the blur kernel points to the correct direction, thereby obtaining a closer approximation The blur kernel of the real scene, its implementation steps include the following:

(1)输入模糊图像y,将模糊图像y设为候选图像;(1) Input a blurred image y, and set the blurred image y as a candidate image;

(2)取大小为3×3的高斯模糊核作为初始化模糊核,用k1表示;(2) Take a Gaussian blur kernel with a size of 3×3 as the initialization blur kernel, denoted by k 1 ;

(3)取全为0的与图像大小相同的二进制掩模作为初始掩模,用M1表示,对外部样例块数据集为BSD500标准数据集进行学习,得到本发明的初始化外部图像样例块;(3) Get the binary mask identical with image size that is all 0 as the initial mask, represent with M 1 , learn for the BSD500 standard data set to the external sample block data set, obtain the initialization external image sample of the present invention piece;

(4)对模糊图像y进行初始化,得到初始候选图像x0 (4) Initialize the blurred image y to get the initial candidate image x 0

xx 00 == argarg mm ii nno xx &Sigma;&Sigma; DD. ** ww ** || || KDKD ** xx -- DD. ** ythe y || || 22 ++ &alpha;&alpha; || || DD. hh xx || || 22 ++ &alpha;&alpha; || || DD. vv xx || || 22 ,, -- -- -- << 11 >>

其中,K代表模糊核k1的矩阵形式,y代表输入的模糊图像,x0代表本次迭代想要得到的清晰候选图像,D*是不同方向上偏微分的矩阵形式,w*是这些不同方向偏微分所对应的标量权重,Dh和Dv分别为水平和垂直方向上的一阶偏导数的矩阵形式,x是和候选图像大小相同的未知矩阵,表示目标函数为最小值时的x的返回值;Among them, K represents the matrix form of the blur kernel k 1 , y represents the input blurred image, x 0 represents the clear candidate image to be obtained in this iteration, D * is the matrix form of partial differential in different directions, and w * is these different The scalar weight corresponding to the direction partial differential, D h and D v are the matrix forms of the first-order partial derivatives in the horizontal and vertical directions, respectively, x is an unknown matrix with the same size as the candidate image, Indicates the return value of x when the objective function is the minimum value;

(5)调用高斯金字塔模型,根据初始化时设定的模糊核k1的大小,计算金字塔总层数N,初始金字塔层数标签t=1;(5) Call the Gaussian pyramid model, according to the size of the fuzzy kernel k 1 set during initialization, calculate the total number of layers N of the pyramid, and the initial pyramid number of layers label t=1;

(6)将候选图像x0根据金字塔层数进行下采样,得到金子塔层第1层的候选图像x1(6) the candidate image x0 is down-sampled according to the number of pyramid layers, and the candidate image x1 of the first layer of the pyramid layer is obtained;

(7)将候选图像xt和模糊核kt根据金子塔层数进行上采样;(7) Upsampling the candidate image x t and the blur kernel k t according to the number of pyramid layers;

(8)判断金字塔标签t是否为N,如果是,保存N层的候选图像xN和模糊核kN执行骤(9),否则执行步骤(14);(8) judge whether the pyramid label t is N, if yes, save the candidate image x N and the fuzzy kernel k N of N layers to perform step (9), otherwise perform step (14);

(9)设置局部迭代最高次数为200,迭代次数标签j=1,将(8)中求得的候选图像xN用Xj表示,作为新的候选图像,将模糊核kN用Kj表示,作为新的模糊核;(9) Set the maximum number of local iterations to 200, the number of iterations label j=1, and represent the candidate image x N obtained in (8) by X j , as a new candidate image, represent the blur kernel k N by K j , as the new blur kernel;

(10)计算当前候选图像Xj的l2范数;(10) Calculate the 12 norm of the current candidate image X j ;

(11)保持模糊核k1以及候选图像Xj的l2范数||Xj||2保持不变,采用l1/l2范数的稀疏正则对图像迭代方向加以限制,根据迭代收缩阈值算法优化公式计算新候选图像Xj+1(11) Keep the blur kernel k 1 and the l 2 norm ||X j || 2 of the candidate image X j unchanged, and use the sparse regularization of the l 1 /l 2 norm to limit the image iteration direction, according to the iterative shrinkage Threshold algorithm optimization formula to calculate new candidate image X j+1 ;

Xx jj ++ 11 == argarg mm ii nno xx &alpha;&alpha; || || xx &CircleTimes;&CircleTimes; KK jj -- ythe y || || 22 22 ++ || || xx || || 11 || || xx || || 22 ++ &beta;&beta; || || KK jj || || 11 ,, -- -- -- << 22 >>

其中,Kj为为j次迭代的模糊核,x为与候选图像大小相同的未知矩阵,y为输入的模糊图像,为二维卷积运算符,式中的第一项为数据保真项,第二项是对x添加的l1/l2范数正则项,最后一项是对模糊核Kj添加的l1范数正则,标量权重α和β用来表示控制模糊核Kj和图像正则项的相对强度,argmin表示目标函数为最小值时的x的值;Among them, K j is the blur kernel of j iterations, x is an unknown matrix with the same size as the candidate image, y is the input blurred image, is a two-dimensional convolution operator, the first item in the formula is the data fidelity item, the second item is the l 1 /l 2 norm regular item added to x, and the last item is l added to the blur kernel K j 1 Norm regularization, scalar weights α and β are used to represent the relative strength of the control blur kernel K j and image regularization items, and argmin represents the value of x when the objective function is the minimum value;

(12)保持候选图像Xj+1不变,根据下式计算新的模糊核Kj+1(12) keep the candidate image X j+1 unchanged, and calculate the new blur kernel K j+1 according to the following formula;

KK jj ++ 11 == argarg mm ii nno kk &alpha;&alpha; || || Xx jj ++ 11 &CircleTimes;&CircleTimes; kk -- ythe y || || 22 22 ++ || || Xx jj ++ 11 || || 11 || || Xx jj ++ 11 || || 22 ,, -- -- -- << 33 >> ``

其中,y为输入的模糊图像,为二维卷积运算符,k为与模糊核大小相同的未知矩阵,arg min表示目标函数为最小值时的k的值,第一项为数据保真项,第二项是对候选图像Xj+1添加的l1/l2范数正则项,标量权重α表示控制模糊核的相对强度,将模糊核求解问题转化为优化问题,采用双共轭梯度解法求解方法,返回函数最小化时的k值,作为新的模糊核Kj+1Among them, y is the input blurred image, is a two-dimensional convolution operator, k is an unknown matrix with the same size as the blur kernel, arg min indicates the value of k when the objective function is the minimum value, the first item is the data fidelity item, and the second item is for the candidate image X The l 1 /l 2 norm regularization item added by j+1 , the scalar weight α represents the relative strength of the control fuzzy kernel, transforms the fuzzy kernel solution problem into an optimization problem, uses the double conjugate gradient solution method to solve the problem, and returns the function when it is minimized The value of k is used as the new blur kernel K j+1 ;

(13)迭代次数标签j加1,重新赋值给j,作为新的迭代次数标签,判断新的迭代次数标签是否为200,如果是,输出候选图像X200以及模糊核K200,否则,返回步骤(10);(13) Add 1 to the iteration number label j, and reassign it to j as a new iteration number label, judge whether the new iteration number label is 200, if yes, output candidate image X 200 and blur kernel K 200 , otherwise, return to step (10);

(14)更新二进制掩模Mt+1:在所有图像块中,计算八个方向的梯度信息,选取边缘信息较强的前2%的图像块,将这些图像块与掩模Mt相对的位置置1,其余位置置0,作为新的二进制掩模Mt+1(14) Update the binary mask M t+1 : in all image blocks, calculate the gradient information in eight directions, select the top 2% image blocks with strong edge information, and compare these image blocks with the mask M t The position is set to 1, and the rest are set to 0, as a new binary mask M t+1 ;

(15)保持二进制掩模Mt+1、外部图像样例块向量Si以及候选图像xt不变,更新图像块的方差ηi(15) Keep the binary mask M t+1 , the external image sample block vector S i and the candidate image x t unchanged, and update the variance η i of the image block;

(16)保持其他参数不变,在二进制掩模Mt+1置1的所有位置,设学习到的图像块为pi=ηiSii,ηi为图像块i的方差,Si为外部图像样例块的向量形式,μi为图像块i的灰度的均值,在外部图像样例块集中找到与候选图像块(Qix-pi)/ηi最相似的样例块Si,得到新的图像样例块Si(16) Keep other parameters unchanged, set the learned image block as p ii S ii at all positions where binary mask M t+1 is set to 1, and η i is the variance of image block i, S i is the vector form of the external image sample block, μ i is the mean value of the gray level of the image block i, find the sample most similar to the candidate image block (Q i xp i )/η i in the external image sample block set block S i , get a new image sample block S i ;

(17)保持其他参数不变,计算得到新的候选图像xt+1(17) keep other parameters unchanged, and calculate a new candidate image x t+1 ;

(18)保持其他参数不变,利用如下公式求解模糊核kt+1 (18) Keeping other parameters unchanged, use the following formula to solve the fuzzy kernel k t+1

kk tt ++ 11 == &Sigma;&Sigma; &delta;&delta; ** &omega;&omega; ** || || kk tt ** &delta;&delta; ** xx tt ++ 11 -- &delta;&delta; ** ythe y || || 22 ++ &beta;&beta; || || kk tt || || 22 -- -- -- << 44 >>

其中δ*代表对应D*的偏导数;y代表输入的模糊图像,w*是这些不同方向偏微分所对应的标量权重,kt表示t金字塔层的模糊核,xt+1是t+1金字塔层的候选图像,设置不在掩模Mt+1中的梯度信息δ*xt为零;Where δ * represents the partial derivative corresponding to D * ; y represents the input blurred image, w * is the scalar weight corresponding to these partial differentials in different directions, k t represents the blur kernel of the t pyramid layer, and x t+1 is t+1 For the candidate image of the pyramid layer, set the gradient information δ * x t not in the mask M t+1 to zero;

(19)金字塔层数标签t加1,重新赋值给t,作为新的金字塔层数标签,返回步骤(7);(19) Pyramid layer number label t adds 1, reassigns to t, as new pyramid number of layer label, returns to step (7);

本发明与现有技术相比具有如下优点:Compared with the prior art, the present invention has the following advantages:

第一,模糊核估计准确First, the blur kernel estimation is accurate

本发明在原有的基于外部图像块先验方法的基础上,在估计模糊核的金字塔模型的最后一层中添加一项模糊核的l1/l2稀疏范数正则项,使得迭代过程中,估计模糊核的阶段指向正确的方向,从而能得到更加逼近真实场景的模糊核,提高了模糊核估计的准确性。On the basis of the original prior method based on external image blocks, the present invention adds an l 1 /l 2 sparse norm regularization term of the fuzzy kernel to the last layer of the pyramid model for estimating the fuzzy kernel, so that in the iterative process, The stage of estimating the blur kernel points to the right direction, so that a blur kernel that is closer to the real scene can be obtained, and the accuracy of blur kernel estimation is improved.

第二,自适应性强Second, strong adaptability

现有的一些图像去模糊的技术,对参数设置要求比较高,参数选择不当时极易出现图像模糊和过拟合的现象,本发明添加l1/l2稀疏范数正则项,不需要过多设置参数,从而使估计模糊核的过程具有较强的自适应性。Some existing image deblurring technologies have relatively high requirements for parameter settings, and image blurring and overfitting are likely to occur when parameters are improperly selected. The present invention adds l 1 /l 2 sparse norm regularization items without over More parameters are set, so that the process of estimating the blur kernel has strong adaptability.

附图说明Description of drawings

图1是本发明的实现流程图;Fig. 1 is the realization flowchart of the present invention;

图2是第一组模糊核估计实验过程图;Fig. 2 is the process diagram of the first group of fuzzy kernel estimation experiments;

图3是第一组实验结果局部放大对比图;Figure 3 is a partial enlarged comparison diagram of the first group of experimental results;

图4是第二组模糊核估计实验过程图;Fig. 4 is a second group of fuzzy kernel estimation experimental process diagrams;

图5是第二组实验结果局部放大对比图;Figure 5 is a partial enlarged comparison diagram of the second group of experimental results;

图6是第三组模糊核估计实验过程图;Fig. 6 is the third group of fuzzy kernel estimation experimental process diagram;

图7是第三组实验结果局部放大对比图;Figure 7 is a partial enlarged comparison diagram of the third group of experimental results;

图8是第四组模糊核估计实验过程图;Fig. 8 is a fourth group of fuzzy kernel estimation experimental process diagrams;

图9是第四组实验结果局部放大对比图。Fig. 9 is a partial enlarged comparison diagram of the fourth group of experimental results.

具体实施方法Specific implementation method

以下参照附图对本发明的技术方案和效果作进一步详细描述。The technical solutions and effects of the present invention will be further described in detail below with reference to the accompanying drawings.

参照图1,本发明的实现步骤如下:With reference to Fig. 1, the realization steps of the present invention are as follows:

步骤1:输入模糊图像y,并将模糊图像y设为候选图像。Step 1: Input a blurred image y, and set the blurred image y as a candidate image.

本实例选取4张各不相同的自然图像,如附图2(a)、4(a)、6(a)、8(a)所示,其名字分别为:brige、Boats、Beverage以及tower,他们的图像尺寸大小分别为:Boats和tower为256×256,图像brige的大小为419×566,图像Beverage的大小为520×395;其中图像Boats为灰度图像,tower、brige以及Beverage为彩色RGB图像。对其进行人工模糊混合处理,得到如附图2(b)、4(b)、6(b)、8(b)所示模糊图像y。In this example, four different natural images are selected, as shown in Figures 2(a), 4(a), 6(a), and 8(a), and their names are: brige, Boats, Beverage, and tower. Their image sizes are: Boats and tower are 256×256, image brige is 419×566, and image Beverage is 520×395; the image Boats is a grayscale image, tower, brige and Beverage are color RGB image. It is artificially blurred and mixed to obtain the blurred image y shown in Figures 2(b), 4(b), 6(b), and 8(b).

步骤2:初始化模糊核Step 2: Initialize the blur kernel

用matlab的fspecial函数生成一个大小为3×3的高斯模糊核作为初始最外层模糊核,用k0表示;Use the fspecial function of matlab to generate a Gaussian blur kernel with a size of 3×3 as the initial outermost blur kernel, denoted by k 0 ;

步骤3:初始化二进制掩模,初始化外部图像样例块。Step 3: Initialize the binary mask and initialize the external image sample block.

取一个大小与图像矩阵相同、数值全为1的矩阵作为初始化二进制掩模,用M1表示;Take a matrix with the same size as the image matrix and all values are 1 as the initial binary mask, denoted by M 1 ;

对外部样例块数据集BSD500标准数据集进行如下学习:The external sample block data set BSD500 standard data set is studied as follows:

首先,对该数据集中所有500张图的每一个维度以采样比例为1:2向下下采样,以初始二进制掩模M1为中心,从500张图像中提取大小为5×5的220K个图像块;First, each dimension of all 500 images in the dataset is down-sampled with a sampling ratio of 1:2, centered on the initial binary mask M 1 , and 220K images of size 5×5 are extracted from the 500 images image blocks;

然后,对这些图像块进行归一化处理,最后设置聚类数目为2560,使用K均值算法对220K图像块进行聚类,形成2560个聚类簇,提取这2560个聚类簇的聚类中心作为需要的外部图像样例块。Then, normalize these image blocks, and finally set the number of clusters to 2560, use the K-means algorithm to cluster 220K image blocks to form 2560 clusters, and extract the cluster centers of these 2560 clusters Sample blocks of external images as needed.

步骤4:对输入模糊图像y进行如下式处理,得到初始化候选图像x0 Step 4: Perform the following processing on the input blurred image y to obtain the initialization candidate image x 0

xx 00 == argarg mm ii nno xx &Sigma;&Sigma; DD. ** ww ** || || KDKD ** xx -- DD. ** ythe y || || 22 ++ &alpha;&alpha; || || DD. hh xx || || 22 ++ &alpha;&alpha; || || DD. vv xx || || 22 ,, -- -- -- << 11 >>

其中,K代表模糊核k1的矩阵形式,y代表输入的模糊图像,x0代表本次迭代想要得到的清晰候选图像,D*是不同方向上偏微分的矩阵形式,w*是这些不同方向偏微分所对应的标量权重,Dh和Dv分别为水平和垂直方向上的一阶偏导数的矩阵形式,x是和候选图像大小相同的未知矩阵,表示目标函数为最小值时的x的返回值||·||2表示矩阵一范数的平方,∑为求和符号。Among them, K represents the matrix form of the blur kernel k 1 , y represents the input blurred image, x 0 represents the clear candidate image to be obtained in this iteration, D * is the matrix form of partial differential in different directions, and w * is these different The scalar weight corresponding to the direction partial differential, D h and D v are the matrix forms of the first-order partial derivatives in the horizontal and vertical directions, respectively, x is an unknown matrix with the same size as the candidate image, Indicates the return value of x when the objective function is the minimum ||||

步骤5:调用高斯金字塔模型,计算金字塔总层数N:Step 5: Call the Gaussian pyramid model to calculate the total number of pyramid layers N:

其中,N为金子塔总层数,为向下取整操作,log表示以2为底的对数操作,b是根据模糊核大小确定的用户参数,初始金字塔层数标签t=1。Among them, N is the total number of layers of the pyramid, It is a rounding down operation, log represents a logarithmic operation with base 2, b is a user parameter determined according to the size of the fuzzy kernel, and the initial pyramid layer number label t=1.

步骤6:对初始化候选图像x0用MATLAB的pyrDown函数进行下采样,将下采样结果作为金子塔层第1层的候选图像x1Step 6: Use MATLAB's pyrDown function to down-sample the initialization candidate image x 0 , and use the down-sampling result as the candidate image x 1 of the first layer of the pyramid.

步骤7:对金子塔第t层的候选图像xt用MATLAB的pyrDown函数进行上采样,并将上采样结果重新赋值给xtStep 7: Use MATLAB's pyrDown function to up-sample the candidate image x t of the t-th layer of the pyramid, and reassign the up-sampling result to x t .

对金子塔第t层的模糊核用MATLAB的pyrDown函数kt进行上采样,并将上采样结果重新赋值给ktUse MATLAB's pyrDown function k t to up-sample the fuzzy kernel of the t-th layer of the pyramid, and reassign the up-sampling result to k t .

步骤8:判断金字塔标签t是否为金字塔层数N,如果是,保存金字塔第N层的候选图像xN和模糊核kN执行骤9,否则执行步骤14;Step 8: Determine whether the pyramid label t is the number of pyramid layers N, if so, save the candidate image x N and the blur kernel k N of the Nth layer of the pyramid and execute step 9, otherwise execute step 14;

步骤9:设置局部迭代最高次数为200,迭代次数标签j=1,将(8)中求得的候选图像xN用Xj表示,作为新的候选图像,模糊核kN用Kj表示,作为新的模糊核。Step 9: Set the maximum number of local iterations to 200, and the number of iterations label j=1. The candidate image x N obtained in (8) is represented by X j as a new candidate image, and the blur kernel k N is represented by K j . as the new blur kernel.

步骤10:计算第j次迭代的候选图像Xj的l2范数||Xj||2Step 10: Calculate the l 2 norm ||X j || 2 of the candidate image X j of the jth iteration;

步骤11:保持模糊核Kj以及候选图像Xj的l2范数||Xj||2保持不变,采用的l1/l2范数的稀疏正则对图像迭代方向加以限制,根据迭代收缩阈值算法优化公式计算新候选图像Xj+1Step 11: Keep the l 2 norm ||X j || 2 of the blur kernel K j and the candidate image X j unchanged, and use the sparse regularization of the l 1 /l 2 norm to limit the image iteration direction, according to the iteration Shrinkage threshold algorithm optimization formula to calculate new candidate image X j+1 ;

Xx jj ++ 11 == argarg mm ii nno xx &alpha;&alpha; || || xx &CircleTimes;&CircleTimes; KK jj -- ythe y || || 22 22 ++ || || xx || || 11 || || xx || || 22 ++ &beta;&beta; || || KK jj || || 11 ,, -- -- -- << 22 >>

其中,Kj为模糊核,x为与候选图像大小相同的未知矩阵,y为输入的模糊图像,为二维卷积运算符。式中的第一项为数据保真项,第二项是对x添加的l1/l2范数正则项,最后一项是对模糊核Kj添加的l1范数正则,α和β为标量权重,用来表示控制模糊核Kj和图像正则项的相对强度,argmin表示目标函数为最小值时的x的值,||·||1表示矩阵一范数,||·|2表示矩阵二范数,表示矩阵二范数的平方;where Kj is the blur kernel, x is an unknown matrix with the same size as the candidate image, y is the input blurred image, is a two-dimensional convolution operator. The first item in the formula is the data fidelity item, the second item is the l 1 /l 2 norm regularization item added to x, and the last item is the l 1 norm regularization addition to the blur kernel K j , α and β is a scalar weight, which is used to represent the relative strength of the control blur kernel K j and the image regularization term, argmin represents the value of x when the objective function is the minimum value, ||·|| 1 represents the matrix-norm, ||·| 2 Represents the second norm of the matrix, Represents the square of the two-norm of the matrix;

步骤12:保持候选图像Xj+1不变,根据下式计算得到第j次迭代的模糊核Kj+1Step 12: keep the candidate image X j+1 unchanged, and calculate the blur kernel K j+1 of the jth iteration according to the following formula;

KK jj ++ 11 == argarg mm ii nno kk &alpha;&alpha; || || Xx jj ++ 11 &CircleTimes;&CircleTimes; kk -- ythe y || || 22 22 ++ || || Xx jj ++ 11 || || 11 || || Xx jj ++ 11 || || 22 ,, -- -- -- << 33 >> ``

其中,Xj+1为候选图像,y为输入的模糊图像,为二维卷积运算符,k为与模糊核大小相同的未知矩阵,arg min表示目标函数为最小值时的k的值;第一项为数据保真项,第二项是对候选图像Xj+1添加的l1/l2范数正则项,标量权重α表示控制模糊核的相对强度,||·||1表示矩阵一范数,||·||2表示矩阵二范数,表示矩阵二范数的平方,将模糊核求解问题转化为优化问题,采用双共轭梯度解法求解方法,返回函数最小化时的k值,作为新的模糊核Kj+1Among them, X j+1 is the candidate image, y is the input blurred image, is a two-dimensional convolution operator, k is an unknown matrix with the same size as the blur kernel, and arg min indicates the value of k when the objective function is the minimum value; the first item is the data fidelity item, and the second item is for the candidate image X The l 1 /l 2 norm regularization item added by j+1 , the scalar weight α represents the relative strength of the control blur kernel, ||·|| 1 represents the matrix one-norm, ||·|| 2 represents the matrix two-norm, Represents the square of the second norm of the matrix, converts the fuzzy kernel solving problem into an optimization problem, adopts the double conjugate gradient solution method to solve the problem, and returns the k value when the function is minimized as a new fuzzy kernel K j+1 ;

步骤13:迭代次数标签j加1,重新赋值给j,作为新的迭代次数标签,判断新的迭代次数标签是否为200,如果是,则输出候选图像X200以及模糊核K200,作为本实例的最终结果,如附图2(c)、4(c)、6(c)、8(c)所示,否则,返回步骤10。Step 13: Add 1 to the iteration number label j, and reassign it to j as a new iteration number label, and judge whether the new iteration number label is 200, and if so, output candidate image X 200 and blur kernel K 200 as this example The final result, as shown in Figures 2(c), 4(c), 6(c), and 8(c), otherwise, return to step 10.

步骤14:更新二进制掩模Mt+1Step 14: Update the binary mask M t+1 .

在所有图像块中,计算八个方向的梯度信息,选取边缘信息较强的前2%的图像块,将这些图像块与掩模Mt相对的位置置1,其余位置置0,作为新的二进制掩模Mt+1In all image blocks, calculate the gradient information in eight directions, select the top 2% image blocks with strong edge information, set the positions of these image blocks relative to the mask M t to 1, and set the rest of the positions to 0, as the new Binary mask M t+1 .

步骤15:保持二进制掩模Mt+1、外部图像样例块向量Si以及候选图像xt不变,用二进制矩阵提取算子Qi提取候选图像x中位置i处大小为5×5个像素的图像块,并更新图像块的方差ηiStep 15: Keep the binary mask M t+1 , the external image sample block vector S i and the candidate image x t unchanged, and use the binary matrix extraction operator Q i to extract the candidate image x at position i with a size of 5×5 The image block of the pixel, and update the variance η i of the image block;

(15a)令vi=Qixt-pi,计算权重系数ωi(15a) Let v i =Q i x t -p i , and calculate the weight coefficient ω i :

&omega;&omega; ii == (( 22 &epsiv;&epsiv; 22 ++ vv ii TT vv ii )) -- 11 -- -- -- << 55 >>

其中,Qi是二进制矩阵提取算子,ε为提前设定的灰度阈值,pi=ηiSii,Si为位置i处对应的外部图像样例块的向量形式,ηi为原图像块的方差值,μi为图像块i的灰度,xt是金字塔第t层的候选图像,是中间变量vi的矩阵转置;Among them, Q i is a binary matrix extraction operator, ε is the gray threshold set in advance, p ii S ii , S i is the vector form of the corresponding external image sample block at position i, η i is the variance value of the original image block, μ i is the gray level of the image block i, x t is the candidate image of the t-th layer of the pyramid, is the matrix transpose of the intermediate variable v i ;

(15b)利用如下公式计算得到新的图像块的方差 (15b) Use the following formula to calculate the variance of the new image block

&eta;&eta; ** ii == &omega;&omega; ii &beta;&beta; || Mm tt ++ 11 || SS ii TT (( QQ ii xx tt -- &mu;&mu; ii )) -- Ff rr ee ff -- 11 (( Ff (( &eta;&eta; ii )) )) &omega;&omega; ii &beta;&beta; || Mm tt ++ 11 || SS ii TT SS ii -- -- -- << 66 >>

其中,Mt+1是二进制掩模,ωi为权重系数,β为正则化图像强度,Si为位置i处对应的外部图像样例块的向量形式,SiT表示Si的矩阵转置,Qi是二进制矩阵提取算子,μi为图像块i的灰度值,|·|表示对矩阵的行列式计算,β为正则化图像强度,xt是金字塔第t层的候选图像,F是当前候选图像方差ηi的经验累积分布,ηi为原图像块的方差值,Fref是外部图像样例块的局部对比参考累计分布。Among them, M t+1 is the binary mask, ω i is the weight coefficient, β is the regularized image intensity, S i is the vector form of the corresponding external image sample block at position i, and S iT represents the matrix transposition of S i , Q i is a binary matrix extraction operator, μ i is the gray value of image block i, |·| represents the determinant calculation of the matrix, β is the regularized image intensity, x t is the candidate image of the t-th layer of the pyramid, F is the empirical cumulative distribution of the variance η i of the current candidate image, η i is the variance value of the original image block, and F ref is the local comparison reference cumulative distribution of the external image sample block.

步骤16:保持其他参数不变,在二进制掩模Mt+1置1的所有位置,设学习到的图像块为pi=ηiSii,ηi为图像块i的方差,Si为外部图像样例块的向量形式,μi为图像块i的灰度的均值,在外部图像样例块集中找到与候选图像块(Qix-pi)/ηi最相似的样例块Si,得到新的图像样例块SiStep 16: keep other parameters unchanged, set the learned image block as p ii S ii at all positions where binary mask M t+1 is set to 1, and η i is the variance of image block i, S i is the vector form of the external image sample block, μ i is the mean value of the gray level of the image block i, find the sample most similar to the candidate image block (Q i xp i )/η i in the external image sample block set block S i , get a new image sample block S i ;

步骤17:保持其他参数不变,计算得到金子塔层第t+1层的候选图像xt+1Step 17: Keeping other parameters unchanged, calculate the candidate image x t+1 of the t+1th layer of the pyramid layer:

xx tt ++ 11 == Ff -- 11 (( BB )) ++ &beta;&beta; || Mm tt ++ 11 || &Sigma;&Sigma; ii &Element;&Element; Mm tt ++ 11 22 22 &epsiv;&epsiv; 22 ++ vv ii TT vv ii QQ ii TT (( &eta;&eta; ii SS ii ++ &mu;&mu; ii )) -- -- -- << 77 >>

B为中间变量,其表示如下:B is an intermediate variable, which is expressed as follows:

式中,xt+1为t+1金字塔层的候选图像,F-1代表傅里叶反变换,β为正则化图像强度,Mt+1是二进制掩模,|·|表示对矩阵的行列式计算,ε为设定的灰度阈值,vi=Qixt-pi,xt是金字塔第t层的候选图像,Qi是一个二进制矩阵提取算子,Qi T表示二进制矩阵提取算子的矩阵转置,ηi为位置i处图像块的方差,μi为图像块的灰度均值,Si为位置i处对应的外部图像样例块的向量形式;在中间变量B中,代表求复共轭运算,δ*代表对应微分矩阵的偏导数⊙为元素级乘法运算符,kt是金字塔t层的模糊核,y为初始输入的模糊图像。In the formula, x t+1 is the candidate image of the t+1 pyramid layer, F -1 represents the inverse Fourier transform, β is the regularized image intensity, M t+1 is the binary mask, and || Determinant calculation, ε is the set gray threshold, v i =Q i x t -p i , x t is the candidate image at the tth layer of the pyramid, Q i is a binary matrix extraction operator, and Q i T represents binary The matrix transposition of the matrix extraction operator, η i is the variance of the image block at position i, μ i is the gray mean value of the image block, S i is the vector form of the corresponding external image sample block at position i; in the intermediate variable B, Represents the complex conjugate operation, δ * represents the partial derivative of the corresponding differential matrix ⊙ is an element-level multiplication operator, k t is the fuzzy kernel of the pyramid layer t, and y is the initial input fuzzy image.

步骤18:保持其他参数不变,利用如下公式求解第t+1金字塔层的模糊核kt+1 Step 18: Keep other parameters unchanged, use the following formula to solve the fuzzy kernel k t+1 of the t+1th pyramid layer

kk tt ++ 11 == &Sigma;&Sigma; &delta;&delta; ** &omega;&omega; ** || || kk tt ** &delta;&delta; ** xx tt ++ 11 -- &delta;&delta; ** ythe y || || 22 ++ &beta;&beta; || || kk tt || || 22 -- -- -- << 44 >>

其中δ*代表对应D*的偏导数,y代表输入的模糊图像,ω*是这些不同方向偏微分所对应的标量权重,kt表示金字塔第t层的模糊核,xt+1是金字塔第t+1层所得的候选图像,设置不在掩模Mt+1中的梯度信息δ*xt+1为零,β表示控制模糊核kt正则项的相对强度,∑为求和符号,||·||表示矩阵一范数,||·||2表示矩阵二范数,表示矩阵二范数的平方。Among them, δ * represents the partial derivative corresponding to D * , y represents the input blurred image, ω * is the scalar weight corresponding to these partial differentials in different directions, k t represents the fuzzy kernel of the t-th layer of the pyramid, and x t+1 is the pyramid's first For the candidate image obtained from layer t+1, set the gradient information δ*x t+1 that is not in the mask M t+1 to zero, β represents the relative strength of the regular term that controls the blur kernel k t , ∑ is the summation symbol, | |·|| represents the matrix one-norm, ||·|| 2 represents the matrix two-norm, Represents the square of the 2-norm of a matrix.

步骤19:给金字塔层数标签t加1,重新赋值给t,作为新的金字塔层数标签,返回步骤7。Step 19: Add 1 to the pyramid layer label t, reassign t as a new pyramid layer label, and return to step 7.

本发明的效果可以通过以下实验来进一步说明:Effect of the present invention can be further illustrated by following experiments:

1、仿真条件:1. Simulation conditions:

本发明的软件环境为Windows 7旗舰版64位系统,MATLAB 2014a。硬件环境为IntelCore2Duo 3.2GHz的CPU,并且内存为DDR34GB的环境下运行。The software environment of the present invention is Windows 7 Ultimate Edition 64-bit system, MATLAB 2014a. The hardware environment is IntelCore2Duo 3.2GHz CPU, and the memory is DDR34GB.

2.仿真内容:2. Simulation content:

仿真1,从已有的非线性模糊核数据集中选取19×19的非线性模糊核对清晰图像图进行模糊混合处理,得到合成的模糊图像图,对合成的模糊图像图进行盲去模糊处理,得到去模糊后的图像和估计出的模糊核图,如图2所示。其中Simulation 1, select 19×19 nonlinear fuzzy kernels from the existing nonlinear fuzzy kernel data set to perform fuzzy mixing processing on the clear image to obtain a synthesized blurred image, and perform blind deblurring on the synthesized blurred image to obtain The deblurred image and the estimated blur kernel map are shown in Figure 2. in

图2(a)为brige的原始清晰图像,图像右下角为合成模糊所用的模糊核;Figure 2(a) is the original clear image of brige, and the lower right corner of the image is the blur kernel used for synthetic blur;

图2(b)为合成的模糊图像;Figure 2(b) is the synthesized blurred image;

图2(c)为去模糊后的图像,图像的右下角为估计出的模糊核;Figure 2(c) is the image after deblurring, and the lower right corner of the image is the estimated blur kernel;

将图2(b)所示的模糊图像与图2(c)去模糊后图像的细节对比,结果如图3。从图3可见,本发明估计出来的模糊核逼近于真实的模糊核,去模糊后图像接近清晰图像。Comparing the blurred image shown in Figure 2(b) with the details of the deblurred image in Figure 2(c), the result is shown in Figure 3. It can be seen from Fig. 3 that the blur kernel estimated by the present invention is close to the real blur kernel, and the image after deblurring is close to the clear image.

仿真2,从已有的非线性模糊核数据集中选取17×17的非线性模糊核对清晰图像图进行模糊混合处理,得到合成的模糊图像图,对合成的模糊图像图进行盲去模糊处理,得到去模糊后的图像和估计出的模糊核图,如图4所示。其中Simulation 2, select a 17×17 nonlinear fuzzy kernel from the existing nonlinear fuzzy kernel data set to perform blurring and mixing processing on the clear image to obtain a synthesized blurred image, and perform blind deblurring on the synthesized blurred image to obtain The deblurred image and the estimated blur kernel map are shown in Figure 4. in

图4(a)为boats的原始清晰图像,图像右下角为合成模糊所用的模糊核;Figure 4(a) is the original clear image of boats, and the lower right corner of the image is the blur kernel used for synthetic blur;

图4(b)为合成的模糊图像;Figure 4(b) is the synthesized blurred image;

图4(c)为去模糊后的图像,图像的右下角为估计出的模糊核;Figure 4(c) is the image after deblurring, and the lower right corner of the image is the estimated blur kernel;

将图4(b)所示的模糊图像与图4(c)去模糊后图像的细节对比,结果如图5。从图5可见,在本发明在该组图中去模糊后图像清晰,估计出来的模糊核几乎与真实的模糊核相同。Comparing the blurred image shown in Figure 4(b) with the details of the deblurred image in Figure 4(c), the result is shown in Figure 5. It can be seen from FIG. 5 that the image is clear after the present invention deblurs the group of images, and the estimated blur kernel is almost the same as the real blur kernel.

仿真3,从已有的非线性模糊核数据集中选取15×15的非线性模糊核对清晰图像图进行模糊混合处理,得到合成的模糊图像图,对合成的模糊图像图进行盲去模糊处理,得到去模糊后的图像和估计出的模糊核图,如图6所示。其中Simulation 3, select a 15×15 nonlinear fuzzy kernel from the existing nonlinear fuzzy kernel data set to perform fuzzy mixing processing on the clear image to obtain a synthesized blurred image, and perform blind deblurring on the synthesized blurred image to obtain The deblurred image and the estimated blur kernel map are shown in Figure 6. in

图6(a)为Beverage的原始清晰图像,图像右下角为合成模糊所用的模糊核;Figure 6(a) is the original clear image of Beverage, and the lower right corner of the image is the blur kernel used for synthetic blur;

图6(b)为合成的模糊图像;Figure 6(b) is the synthesized blurred image;

图6(c)为去模糊后的图像,图像的右下角为估计出的模糊核;Figure 6(c) is the image after deblurring, and the lower right corner of the image is the estimated blur kernel;

将图6(b)所示的模糊图像与图6(c)去模糊后图像的细节对比,结果如图7。从图7可见,本发明去模糊后图像清晰,估计出来的模糊核几乎与真实的模糊核相同。Comparing the blurred image shown in Figure 6(b) with the details of the deblurred image in Figure 6(c), the result is shown in Figure 7. It can be seen from Fig. 7 that the image is clear after deblurring by the present invention, and the estimated blur kernel is almost the same as the real blur kernel.

仿真4,从已有的非线性模糊核数据集中选取13×13的非线性模糊核对清晰图像图进行模糊混合处理,得到合成的模糊图像图,对合成的模糊图像图进行盲去模糊处理,得到去模糊后的图像和估计出的模糊核图,如图8所示。其中Simulation 4. Select 13×13 nonlinear fuzzy kernels from the existing nonlinear fuzzy kernel data set to perform fuzzy mixing processing on the clear image to obtain a synthesized blurred image, and perform blind deblurring processing on the synthesized blurred image to obtain The deblurred image and the estimated blur kernel map are shown in Figure 8. in

图8(a)为Beverage的原始清晰图像,图像右下角为合成模糊所用的模糊核;Figure 8(a) is the original clear image of Beverage, and the lower right corner of the image is the blur kernel used for synthetic blur;

图8(b)为合成的模糊图像;Figure 8(b) is the synthesized blurred image;

图8(c)为去模糊后的图像,图像的右下角为估计出的模糊核;Figure 8(c) is the image after deblurring, and the lower right corner of the image is the estimated blur kernel;

将图8(b)所示的模糊图像与图8(c)去模糊后图像的细节对比,结果如图9。从图9可见,本发明去模糊后图像清晰,估计出来的模糊核几乎与真实的模糊核相同。Comparing the blurred image shown in Figure 8(b) with the details of the deblurred image in Figure 8(c), the result is shown in Figure 9. It can be seen from Fig. 9 that the image is clear after deblurring by the present invention, and the estimated blur kernel is almost the same as the real blur kernel.

3、实验结果分析:3. Analysis of experimental results:

综合四组实验结果可以看出,本发明去模糊后图像清晰,估计出的模糊核逼近于真实模糊的模糊核,并且一定程度上克服了模糊核中存在类似噪声的现象,体现了本发明对模糊核的估计具有很好的准确性。四组实验的模糊核估计结果表现稳定,体现了本发明具有很好的适应性和稳定性。Comprehensive four groups of experimental results can be seen, the image after deblurring of the present invention is clear, and the fuzzy kernel estimated is close to the fuzzy kernel of real blur, and overcomes to a certain extent the phenomenon that there is similar noise in the fuzzy kernel, has reflected the present invention to The estimation of the blur kernel has good accuracy. The fuzzy kernel estimation results of the four groups of experiments are stable, which shows that the present invention has good adaptability and stability.

Claims (3)

1.一种基于图像块先验与稀疏范数的盲去模糊方法,包括: 1. A blind deblurring method based on image block prior and sparse norm, comprising: (1)输入模糊图像y,将模糊图像y设为候选图像; (1) Input a blurred image y, and set the blurred image y as a candidate image; (2)取大小为3×3的高斯模糊核作为初始化模糊核,用k1表示; (2) Take a Gaussian blur kernel with a size of 3×3 as the initialization blur kernel, denoted by k 1 ; (3)取全为0的与图像大小相同的二进制掩模作为初始掩模,用M1表示,对外部样例块数据集为BSD500标准数据集进行学习,得到本发明的初始化外部图像样例块; (3) Get the binary mask identical with image size that is all 0 as the initial mask, represent with M 1 , learn for the BSD500 standard data set to the external sample block data set, obtain the initialization external image sample of the present invention piece; (4)对模糊图像y进行初始化,得到初始候选图像x0 (4) Initialize the blurred image y to get the initial candidate image x 0 其中,k1代表模糊核k1的矩阵形式,y代表输入的模糊图像,x0代表本次迭代想要得到的清晰候选图像,D*是不同方向上偏微分的矩阵形式,w*是这些不同方向偏微分所对应的标量权重,Dh和Dv分别为水平和垂直方向上的一阶偏导数的矩阵形式,x是和候选图像大小相同的未知矩阵,表示目标函数为最小值时的x的返回值; Among them, k 1 represents the matrix form of the blur kernel k 1 , y represents the input blurred image, x 0 represents the clear candidate image that this iteration wants to obtain, D * is the matrix form of partial differential in different directions, and w * is these The scalar weights corresponding to the partial differentials in different directions, D h and D v are the matrix forms of the first-order partial derivatives in the horizontal and vertical directions, respectively, x is an unknown matrix with the same size as the candidate image, Indicates the return value of x when the objective function is the minimum value; (5)调用高斯金字塔模型,根据初始化时设定的模糊核k1的大小,计算金字塔总层数N,初始金字塔层数标签t=1; (5) Call the Gaussian pyramid model, according to the size of the fuzzy kernel k 1 set during initialization, calculate the total number of layers N of the pyramid, and the initial pyramid number of layers label t=1; (6)将候选图像x0根据金字塔层数进行下采样,得到金子塔层第1层的候选图像x1(6) the candidate image x0 is down-sampled according to the number of pyramid layers, and the candidate image x1 of the first layer of the pyramid layer is obtained; (7)将候选图像xt和模糊核kt根据金子塔层数进行上采样; (7) Upsampling the candidate image x t and the blur kernel k t according to the number of pyramid layers; (8)判断金字塔标签t是否为N,如果是,保存N层的候选图像xN和模糊核kN执行骤(9),否则执行步骤(14); (8) judge whether the pyramid label t is N, if yes, save the candidate image x N and the fuzzy kernel k N of N layers to perform step (9), otherwise perform step (14); (9)设置局部迭代最高次数为200,迭代次数标签j=1,将(8)中求得的候选图像xN用Xj表示,作为新的候选图像,将模糊核kN用Kj表示,作为新的模糊核; (9) Set the maximum number of local iterations to 200, the number of iterations label j=1, and represent the candidate image x N obtained in (8) by X j , as a new candidate image, represent the blur kernel k N by K j , as the new blur kernel; (10)计算当前候选图像Xj的l2范数; (10) Calculate the 12 norm of the current candidate image X j ; (11)保持模糊核k1以及候选图像Xj的l2范数||Xj||2保持不变,采用l1/l2范数的稀疏正则对图像迭代方向加以限制,根据迭代收缩阈值算法优化公式计算新候选图像Xj+1(11) Keep the blur kernel k 1 and the l 2 norm ||X j || 2 of the candidate image X j unchanged, and use the sparse regularization of the l 1 /l 2 norm to limit the image iteration direction, according to the iterative shrinkage Threshold algorithm optimization formula to calculate new candidate image X j+1 ; 其中,Kj为为j次迭代的模糊核,x为与候选图像大小相同的未知矩阵,y为输入的模糊图像,为二维卷积运算符,式中的第一项为数据保真项,第二项是对x添加的l1/l2范数正则项,最后一项是对模糊核Kj添加的l1范数正则,标量权重α和β用来表示控制模糊核Kj和图像正则项的相对强度,argmin表示目标函数为最小值时的x的值; Among them, K j is the blur kernel of j iterations, x is an unknown matrix with the same size as the candidate image, y is the input blurred image, is a two-dimensional convolution operator, the first item in the formula is the data fidelity item, the second item is the l 1 /l 2 norm regular item added to x, and the last item is l added to the blur kernel K j 1 Norm regularization, scalar weights α and β are used to represent the relative strength of the control blur kernel K j and image regularization items, and argmin represents the value of x when the objective function is the minimum value; (12)保持候选图像Xj+1不变,根据下式计算新的模糊核Kj+1(12) keep the candidate image X j+1 unchanged, and calculate the new blur kernel K j+1 according to the following formula; 其中,y为输入的模糊图像,为二维卷积运算符,k为与模糊核大小相同的未知矩阵,arg min表示目标函数为最小值时的k的值,第一项为数据保真项,第二项是对候选图像Xj+1添加的l1/l2范数正则项,标量权重α表示控制模糊核的相对强度,将模糊核求解问题转化为优化问题,采用双共轭梯度解法求解方法,返回函数最小化时的k值,作为新的模糊核Kj+1Among them, y is the input blurred image, is a two-dimensional convolution operator, k is an unknown matrix with the same size as the blur kernel, arg min indicates the value of k when the objective function is the minimum value, the first item is the data fidelity item, and the second item is for the candidate image X The l 1 /l 2 norm regularization item added by j+1 , the scalar weight α represents the relative strength of the control fuzzy kernel, transforms the fuzzy kernel solution problem into an optimization problem, uses the double conjugate gradient solution method to solve the problem, and returns the function when it is minimized The value of k is used as the new blur kernel K j+1 ; (13)迭代次数标签j加1,重新赋值给j,作为新的迭代次数标签,判断新的迭代次数标签是否为200,如果是,输出候选图像X200以及模糊核K200,否则,返回步骤(10); (13) Add 1 to the iteration number label j, and reassign it to j as a new iteration number label, judge whether the new iteration number label is 200, if yes, output candidate image X 200 and blur kernel K 200 , otherwise, return to step (10); (14)更新二进制掩模Mt+1:在所有图像块中,计算八个方向的梯度信息,选取边缘信息较强的前2%的图像块,将这些图像块与掩模Mt相对的位置置1,其余位置置0,作为新的二进制掩模Mt+1(14) Update the binary mask M t+1 : in all image blocks, calculate the gradient information in eight directions, select the top 2% image blocks with strong edge information, and compare these image blocks with the mask M t The position is set to 1, and the rest are set to 0, as a new binary mask M t+1 ; (15)保持二进制掩模Mt+1、外部图像样例块向量Si以及候选图像xt不变,更新图像块的方差ηi(15) Keep the binary mask M t+1 , the external image sample block vector S i and the candidate image x t unchanged, and update the variance η i of the image block; (16)保持其他参数不变,在二进制掩模Mt+1置1的所有位置,设学习到的图像块为pi=ηiSii,ηi为图像块i的方差,Si为外部图像样例块的向量形式,μi为图像块i的灰度的均值,在外部图像样例块集中找到与候选图像块(Qix-pi)/ηi最相似的样例块Si,得到新的图像样例块Si(16) Keep other parameters unchanged, set the learned image block as p ii S ii at all positions where binary mask M t+1 is set to 1, and η i is the variance of image block i, S i is the vector form of the external image sample block, μ i is the mean value of the gray level of the image block i, find the sample most similar to the candidate image block (Q i xp i )/η i in the external image sample block set block S i , get a new image sample block S i ; (17)保持其他参数不变,计算得到新的候选图像xt+1(17) keep other parameters unchanged, and calculate a new candidate image x t+1 ; (18)保持其他参数不变,利用如下公式求解模糊核kt+1 (18) Keeping other parameters unchanged, use the following formula to solve the fuzzy kernel k t+1 其中δ*代表对应D*的偏导数;y代表输入的模糊图像,w*是这些不同方向偏微分所对应的标量权重,kt表示t金字塔层的模糊核,xt+1是t+1金字塔层的候选图像,设置不在掩模Mt+1中的梯度信息δ*xt为零; Where δ * represents the partial derivative corresponding to D * ; y represents the input blurred image, w * is the scalar weight corresponding to these partial differentials in different directions, k t represents the blur kernel of the t pyramid layer, and x t+1 is t+1 For the candidate image of the pyramid layer, set the gradient information δ * x t not in the mask M t+1 to zero; (19)金字塔层数标签t加1,重新赋值给t,作为新的金字塔层数标签,返回步骤(7)。 (19) Add 1 to the pyramid layer label t, and reassign it to t as a new pyramid layer label, and return to step (7). 2.根据权利要求1所述的基于图像块先验与稀疏范数的盲去模糊方法,其中步骤(15)中计算得到新图像块的方差按如下步骤进行: 2. the blind deblurring method based on image block prior and sparse norm according to claim 1, wherein the variance of new image block is calculated in step (15) Proceed as follows: (15a)令vi=Qixt-pi,计算权重系数ωi(15a) Let v i =Q i x t -p i , and calculate the weight coefficient ω i : 其中,Qi是二进制矩阵提取算子,ε为提前设定的灰度阈值,pi=ηiSii,Si为位置i处对应的外部图像样例块的向量形式,ηi为原图像块的方差值,μi为图像块i的灰度值; Among them, Q i is a binary matrix extraction operator, ε is the gray threshold value set in advance, p ii S ii , S i is the vector form of the corresponding external image sample block at position i, η i is the variance value of the original image block, μ i is the gray value of the image block i; (15b)利用如下公式计算得到新的图像块的方差η* i(15b) Utilize the following formula to calculate the variance η * i of the new image block: 其中,Mt+1是二进制掩模,ηi为原图像块的方差值,Qi是二进制矩阵提取算子,Si为位置i处对应的外部图像样例块的向量形式,μi为图像块i的灰度值,β为正则化图像强度,xt是金字塔第t层的候选图像,F是当前候选图像方差ηi的经验累积分布,Fref是外部图像样例块的局部对比参考累计分布。 Among them, M t+1 is a binary mask, η i is the variance value of the original image block, Q i is a binary matrix extraction operator, S i is the vector form of the external image sample block corresponding to position i, μ i is the gray value of the image block i, β is the regularized image intensity, x t is the candidate image at the tth layer of the pyramid, F is the empirical cumulative distribution of the variance η i of the current candidate image, F ref is the local sample block of the external image Compare the reference cumulative distribution. 3.根据权利要求1所述的基于图像块先验与稀疏范数的盲去模糊方法,其中步骤(17)计算新的候选图像xt+1,通过如下公式进行: 3. The blind deblurring method based on image block prior and sparse norm according to claim 1, wherein step (17) calculates a new candidate image x t+1 by the following formula: 其中,xt+1为t+1金字塔层的候选图像,F-1代表傅里叶反变换,β为正则化图像强度,Mt+1是二进制掩模,ε为设定的灰度阈值,vi=Qixt-pi,xt是金字塔第t层的候选图像,Qi是一个二进制矩阵提取算子,Qi T表示二进制矩阵提取算子的矩阵转置,ηi为位置i处图像块的方差,μi为图像块的灰度均值,Si为位置i处对应的外部图像样例块的向量形式;在中间变量B中,代表求复共轭运算,δ*代表对应微分矩阵的偏导数,⊙为元素级乘法运算符,kt是金字塔t层的模糊核,y为初始输入的模糊图像。 Among them, x t+1 is the candidate image of the t+1 pyramid layer, F -1 represents the inverse Fourier transform, β is the regularized image intensity, M t+1 is the binary mask, and ε is the set gray threshold , v i =Q i x t -p i , x t is the candidate image of the tth layer of the pyramid, Q i is a binary matrix extraction operator, Q i T represents the matrix transposition of the binary matrix extraction operator, and η i is The variance of the image block at position i, μ i is the gray mean value of the image block, S i is the vector form of the corresponding external image sample block at position i; in the intermediate variable B, Represents the complex conjugate operation, δ * represents the partial derivative of the corresponding differential matrix, ⊙ is the element-level multiplication operator, k t is the blur kernel of the pyramid layer t, and y is the initial input blur image.
CN201610248012.5A 2016-04-20 2016-04-20 Blind deblurring method based on image block prior and sparse norm Active CN105957024B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201610248012.5A CN105957024B (en) 2016-04-20 2016-04-20 Blind deblurring method based on image block prior and sparse norm

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201610248012.5A CN105957024B (en) 2016-04-20 2016-04-20 Blind deblurring method based on image block prior and sparse norm

Publications (2)

Publication Number Publication Date
CN105957024A true CN105957024A (en) 2016-09-21
CN105957024B CN105957024B (en) 2019-06-18

Family

ID=56917818

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201610248012.5A Active CN105957024B (en) 2016-04-20 2016-04-20 Blind deblurring method based on image block prior and sparse norm

Country Status (1)

Country Link
CN (1) CN105957024B (en)

Cited By (13)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106373107A (en) * 2016-12-06 2017-02-01 四川长虹电器股份有限公司 Automatic image deblurring system and automatic image deblurring method of smart phone
CN106530261A (en) * 2016-12-28 2017-03-22 同观科技(深圳)有限公司 Double-dynamic blurred image restoration method
CN106934216A (en) * 2017-02-16 2017-07-07 山东大学齐鲁医院 Medicine equipment clinical evaluation method based on multiple target
CN107292836A (en) * 2017-06-02 2017-10-24 河海大学常州校区 Image Blind deblurring method based on external image block prior information and rarefaction representation
CN107680062A (en) * 2017-10-12 2018-02-09 长沙全度影像科技有限公司 A kind of micro- burnt Restoration method of blurred image based on l1/l2 priori combination Gaussian priors
CN110084762A (en) * 2019-04-26 2019-08-02 华南理工大学 A kind of deep learning is against convolution model outlier processing method
WO2019174068A1 (en) * 2018-03-15 2019-09-19 华中科技大学 Distance-weighted sparse representation priori-based image restoration and matching integration method
CN110490822A (en) * 2019-08-11 2019-11-22 浙江大学 Method and device for image motion blur removal
CN110728626A (en) * 2018-07-16 2020-01-24 宁波舜宇光电信息有限公司 Image deblurring method and apparatus and training thereof
CN113793272A (en) * 2021-08-11 2021-12-14 东软医疗系统股份有限公司 Image noise reduction method and device, storage medium and terminal
CN114359236A (en) * 2022-01-07 2022-04-15 北京微纳星空科技有限公司 Fuzzy kernel determining method, device, equipment and storage medium of image processing model
CN114998100A (en) * 2020-06-08 2022-09-02 广州超视计生物科技有限公司 Image processing system and method
CN116091367A (en) * 2023-04-10 2023-05-09 中国科学院空天信息创新研究院 Blind deblurring method, device, equipment and medium for optical remote sensing image

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104599242A (en) * 2014-12-09 2015-05-06 西安电子科技大学 Multi-scale non-local regularization blurring kernel estimation method
KR101544171B1 (en) * 2014-12-29 2015-08-12 연세대학교 산학협력단 Apparatus and Method for Super Resolution using Hybrid Feature
CN105184744A (en) * 2015-08-24 2015-12-23 西安电子科技大学 Prior fuzzy kernel estimation method based on standardized sparse measurement image block

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104599242A (en) * 2014-12-09 2015-05-06 西安电子科技大学 Multi-scale non-local regularization blurring kernel estimation method
KR101544171B1 (en) * 2014-12-29 2015-08-12 연세대학교 산학협력단 Apparatus and Method for Super Resolution using Hybrid Feature
CN105184744A (en) * 2015-08-24 2015-12-23 西安电子科技大学 Prior fuzzy kernel estimation method based on standardized sparse measurement image block

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
WEISHENG DONG 等: "Image Deblurring and Super-Resolution by Adaptive Sparse Domain Selection and Adaptive Regularization", 《IEEE TRANSACTIONS ON IMAGE PROCESSING》 *
蔡德生 等: "基于字典稀疏表示和梯度稀疏的图像盲去模糊", 《燕山大学学报》 *

Cited By (20)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106373107A (en) * 2016-12-06 2017-02-01 四川长虹电器股份有限公司 Automatic image deblurring system and automatic image deblurring method of smart phone
CN106373107B (en) * 2016-12-06 2019-03-19 四川长虹电器股份有限公司 Smart phone automated graphics deblurring system and method
CN106530261A (en) * 2016-12-28 2017-03-22 同观科技(深圳)有限公司 Double-dynamic blurred image restoration method
CN106530261B (en) * 2016-12-28 2019-03-19 同观科技(深圳)有限公司 A kind of double dynamic fuzzy image recovery methods
CN106934216A (en) * 2017-02-16 2017-07-07 山东大学齐鲁医院 Medicine equipment clinical evaluation method based on multiple target
CN107292836A (en) * 2017-06-02 2017-10-24 河海大学常州校区 Image Blind deblurring method based on external image block prior information and rarefaction representation
CN107292836B (en) * 2017-06-02 2020-06-26 河海大学常州校区 A Blind Image Deblurring Method Based on External Image Patch Prior Information and Sparse Representation
CN107680062A (en) * 2017-10-12 2018-02-09 长沙全度影像科技有限公司 A kind of micro- burnt Restoration method of blurred image based on l1/l2 priori combination Gaussian priors
WO2019174068A1 (en) * 2018-03-15 2019-09-19 华中科技大学 Distance-weighted sparse representation priori-based image restoration and matching integration method
CN110728626A (en) * 2018-07-16 2020-01-24 宁波舜宇光电信息有限公司 Image deblurring method and apparatus and training thereof
CN110084762A (en) * 2019-04-26 2019-08-02 华南理工大学 A kind of deep learning is against convolution model outlier processing method
CN110084762B (en) * 2019-04-26 2022-11-18 华南理工大学 Deep learning inverse convolution model abnormal value processing method
CN110490822A (en) * 2019-08-11 2019-11-22 浙江大学 Method and device for image motion blur removal
CN110490822B (en) * 2019-08-11 2022-02-15 浙江大学 Method and device for removing motion blur of image
CN114998100A (en) * 2020-06-08 2022-09-02 广州超视计生物科技有限公司 Image processing system and method
CN113793272A (en) * 2021-08-11 2021-12-14 东软医疗系统股份有限公司 Image noise reduction method and device, storage medium and terminal
CN113793272B (en) * 2021-08-11 2024-01-26 东软医疗系统股份有限公司 Image noise reduction method and device, storage medium and terminal
CN114359236A (en) * 2022-01-07 2022-04-15 北京微纳星空科技有限公司 Fuzzy kernel determining method, device, equipment and storage medium of image processing model
CN116091367A (en) * 2023-04-10 2023-05-09 中国科学院空天信息创新研究院 Blind deblurring method, device, equipment and medium for optical remote sensing image
CN116091367B (en) * 2023-04-10 2023-07-18 中国科学院空天信息创新研究院 Optical remote sensing image blind deblurring method, device, equipment and medium

Also Published As

Publication number Publication date
CN105957024B (en) 2019-06-18

Similar Documents

Publication Publication Date Title
CN105957024A (en) Blind deblurring method based on image block prior and sparse norm
Köhler et al. Robust multiframe super-resolution employing iteratively re-weighted minimization
US9042672B2 (en) Motion deblurring using image upsampling
CN102326379B (en) Method for removing blur from image
Wu et al. Learning-based super resolution using kernel partial least squares
CN111369487B (en) Hyperspectral and multispectral image fusion method, system and medium
CN103150713B (en) Utilize the image super-resolution method that image block classification rarefaction representation is polymerized with self-adaptation
US20170109873A1 (en) Image enhancement using self-examples and external examples
Wang et al. Image deblurring with filters learned by extreme learning machine
CN106251297A (en) A kind of estimation based on multiple image fuzzy core the rebuilding blind super-resolution algorithm of improvement
Gao et al. Bayesian image super-resolution with deep modeling of image statistics
Gao et al. Stacked convolutional auto-encoders for single space target image blind deconvolution
Li et al. Single image super-resolution using combined total variation regularization by split Bregman Iteration
Li et al. A maximum a posteriori estimation framework for robust high dynamic range video synthesis
Huang et al. Fast blind image super resolution using matrix-variable optimization
Pires et al. Image denoising using attention-residual convolutional neural networks
Zhang et al. Fast restoration of star image under dynamic conditions via lp regularized intensity prior
Ghassab et al. Plug-and-play video super-resolution using edge-preserving filtering
Li et al. A mixed noise removal algorithm based on multi-fidelity modeling with nonsmooth and nonconvex regularization
CN106033595B (en) A Blind Image Deblurring Method Based on Local Constraints
Lu et al. Smoothing denoising method of spatial filtering image based on Contourlet transform
Dasari et al. Image Enhancement of Underwater Images Using Deep Learning Techniques
CN114418877B (en) Image non-blind deblurring method and system based on gradient amplitude similarity
CN109948571A (en) An optical remote sensing image ship detection method
Xian et al. Hybrid example-based single image super-resolution

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant