Engee 文档
Notebook

模糊图像恢复示例中具有约束的线性最小二乘法

导言

在这个例子中,我们来看看如何通过使用带约束的线性最小二乘法求解优化问题来恢复模糊图像。

我们将附加必要的库。

In [ ]:
# import Pkg
# Pkg.add(["FileIO", "Images", "ImageShow", "LinearAlgebra", "SparseArrays", "JuMP", "Ipopt", "Printf"])
using FileIO, Images, ImageShow, LinearAlgebra, SparseArrays, JuMP, Ipopt, Printf

设置任务

有一张照片显示了Engee展台。

In [ ]:
P = load("image.jpg")
P = Gray.(P)
P = Float64.(channelview(P))
m, n = size(P)
mn = m * n
p1 = heatmap(reverse(P, dims=1), aspect_ratio=:equal, color=:grays, axis=false, colorbar=false, title="原始图像是$m x$n($mn像素)")
display(p1)

任务是拍摄这张照片的模糊版本,并尝试恢复其清晰度。 原始图像是黑白的,这意味着它可以表示为矩阵。 大小 ,由像素值范围从 以前 .

模糊图像

让我们通过在顶部和底部用15个像素平均每个像素来模拟垂直运动模糊的效果。 让我们构造一个稀疏矩阵D,通过矩阵乘法来执行模糊。

In [ ]:
blur = 15
weight = 1/(2*blur + 1)
mindex = Int[]
nindex = Int[]
values = Float64[]
for i = 1:mn
    push!(mindex, i)
    push!(nindex, i)
    push!(values, weight)
    for k = 1:blur
        if i + k <= mn
            push!(mindex, i + k)
            push!(nindex, i)
            push!(values, weight)
        end
        
        if i - k >= 1
            push!(mindex, i - k)
            push!(nindex, i)
            push!(values, weight)
        end
    end
end

D = sparse(mindex, nindex, values, mn, mn);

让我们显示矩阵D的前30×30个元素。

In [ ]:
println("\Blur迭代器(前30x30个元素):")
for i in 1:30
    for j in 1:30
        if abs(i-j) <= blur
            print("x ")
        else
            print("0 ")
        end
    end
    println("...")
end
println("..."^10)
Оператор размытия (первые 30x30 элементов):
x x x x x x x x x x x x x x x x 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ...
x x x x x x x x x x x x x x x x x 0 0 0 0 0 0 0 0 0 0 0 0 0 ...
x x x x x x x x x x x x x x x x x x 0 0 0 0 0 0 0 0 0 0 0 0 ...
x x x x x x x x x x x x x x x x x x x 0 0 0 0 0 0 0 0 0 0 0 ...
x x x x x x x x x x x x x x x x x x x x 0 0 0 0 0 0 0 0 0 0 ...
x x x x x x x x x x x x x x x x x x x x x 0 0 0 0 0 0 0 0 0 ...
x x x x x x x x x x x x x x x x x x x x x x 0 0 0 0 0 0 0 0 ...
x x x x x x x x x x x x x x x x x x x x x x x 0 0 0 0 0 0 0 ...
x x x x x x x x x x x x x x x x x x x x x x x x 0 0 0 0 0 0 ...
x x x x x x x x x x x x x x x x x x x x x x x x x 0 0 0 0 0 ...
x x x x x x x x x x x x x x x x x x x x x x x x x x 0 0 0 0 ...
x x x x x x x x x x x x x x x x x x x x x x x x x x x 0 0 0 ...
x x x x x x x x x x x x x x x x x x x x x x x x x x x x 0 0 ...
x x x x x x x x x x x x x x x x x x x x x x x x x x x x x 0 ...
x x x x x x x x x x x x x x x x x x x x x x x x x x x x x x ...
x x x x x x x x x x x x x x x x x x x x x x x x x x x x x x ...
0 x x x x x x x x x x x x x x x x x x x x x x x x x x x x x ...
0 0 x x x x x x x x x x x x x x x x x x x x x x x x x x x x ...
0 0 0 x x x x x x x x x x x x x x x x x x x x x x x x x x x ...
0 0 0 0 x x x x x x x x x x x x x x x x x x x x x x x x x x ...
0 0 0 0 0 x x x x x x x x x x x x x x x x x x x x x x x x x ...
0 0 0 0 0 0 x x x x x x x x x x x x x x x x x x x x x x x x ...
0 0 0 0 0 0 0 x x x x x x x x x x x x x x x x x x x x x x x ...
0 0 0 0 0 0 0 0 x x x x x x x x x x x x x x x x x x x x x x ...
0 0 0 0 0 0 0 0 0 x x x x x x x x x x x x x x x x x x x x x ...
0 0 0 0 0 0 0 0 0 0 x x x x x x x x x x x x x x x x x x x x ...
0 0 0 0 0 0 0 0 0 0 0 x x x x x x x x x x x x x x x x x x x ...
0 0 0 0 0 0 0 0 0 0 0 0 x x x x x x x x x x x x x x x x x x ...
0 0 0 0 0 0 0 0 0 0 0 0 0 x x x x x x x x x x x x x x x x x ...
0 0 0 0 0 0 0 0 0 0 0 0 0 0 x x x x x x x x x x x x x x x x ...
..............................

将矩阵D乘以变换为列向量的原始图像,得到模糊图像。

In [ ]:
G_vec = D * vec(P)
G = reshape(G_vec, m, n)
p2 = heatmap(reverse(G, dims=1), aspect_ratio=:equal, color=:grays, axis=false, colorbar=false, title="模糊图像")
display(p2)

图像变得模糊,标签不可读。

图像恢复

为了恢复清晰度,假设模糊运算符是已知的。 . 模糊消除和原始图像恢复的效果如何?

最简单的方法是解决最小二乘问题 :

有条件的 .

在这个问题中,模糊矩阵 被设置,并且需要找到这样的向量 ,其中 它将尽可能接近 . 为了使解表示有意义的像素值,将其限制在从 以前 .

In [ ]:
model = Model(Ipopt.Optimizer)
set_silent(model)
@variable(model, 0 <= x[1:mn] <= 1)
@expression(model, expr, D * x - G_vec)
@objective(model, Min, expr' * expr)

我们来解决优化问题。

In [ ]:
optimize!(model)
x_val = value.(x)
xpic = reshape(x_val, m, n)

我们将显示恢复的图像。

In [ ]:
println("\恢复的图像")
p3 = heatmap(reverse(xpic, dims=1), aspect_ratio=:equal, color=:grays, axis=false, colorbar=false, title="恢复的图像")
display(p3)
Восстановленное изображение

通过应用优化问题解决方案,我们能够几乎完全恢复原始图像,使其上的铭文再次变得可区分和可读。

结论

本文以消除垂直模糊为例,考虑了图像恢复的实际任务。 实验结果令人信服地证明了线性最小二乘法在图像处理中求解逆问题的有效性,清楚地表明抽象的数学概念是如何在解决数字视觉信息处理的实际问题中得到具体应用的。