diff --git a/scripts/run.jl b/scripts/run.jl
index e7a2ed44c0575a68a380081e867b9f6e86866fc2..97bb1781224efc526324bf1ee998b07e16a91009 100644
--- a/scripts/run.jl
+++ b/scripts/run.jl
@@ -27,3 +27,5 @@ ctx(experiment_denoise, "fem/denoise")
ctx(experiment_inpaint, "fem/inpaint")
# adaptive optical flow
ctx(experiment_optflow_middlebury_all, "fem/optflow/middlebury")
+
+# additional experiments for paper
diff --git a/scripts/run_experiments.jl b/scripts/run_experiments.jl
index 104a5e861b239dfe24a8dbaf187b9fa1e42b4d16..94670c932dfe694e39a6797d3717273f2fe7d6c8 100644
--- a/scripts/run_experiments.jl
+++ b/scripts/run_experiments.jl
@@ -28,6 +28,38 @@ isdefined(Main, :Revise) && Revise.track(joinpath(@__DIR__, "util.jl"))
using .Util
+module Operators
+
+struct TIdentity end
+(::TIdentity)(tdata, u) = u
+(::TIdentity)(::typeof(adjoint), tdata, v) = v
+Base.show(io, ::MIME"text/latex", ::TIdentity) =
+ print(io, raw"\ensuremath{I}")
+
+struct TInpaint end
+(::TInpaint)(tdata, u) = abs(tdata[begin] - 1.) < 1e-8 ? u : zero(u)
+(T::TInpaint)(::typeof(adjoint), tdata, v) = T(tdata, v)
+Base.show(io, ::MIME"text/latex", ::TInpaint) =
+ print(io, raw"\ensuremath{T_{\mathrm{inpaint}}}")
+
+struct TOptflow end
+(::TOptflow)(tdata, u) = tdata * u
+(::TOptflow)(::typeof(adjoint), tdata, v) = tdata' * v
+Base.show(io, ::MIME"text/latex", ::TOptflow) =
+ print(io, raw"\ensuremath{T_{\mathrm{optflow}}}")
+
+struct SIdentity end
+(::SIdentity)(u, nablau) = u
+Base.show(io, ::MIME"text/latex", ::SIdentity) =
+ print(io, raw"\ensuremath{I}")
+
+struct SNabla end
+(::SNabla)(u, nablau) = nablau'
+Base.show(io, ::MIME"text/latex", ::SNabla) =
+ print(io, raw"\ensuremath{\nabla}")
+
+end
+
grayclamp(value) = Gray(clamp(value, 0., 1.))
loadimg(x) = Float64.(FileIO.load(x))
saveimg(io, x::Array{<:Gray}) = FileIO.save(io, grayclamp.(x))
@@ -156,7 +188,7 @@ function OptFlowState(mesh;
# tdata will be something like nabla(fw)
T(tdata, u) = tdata * u
T(::typeof(adjoint), tdata, v) = tdata' * v
- S(u, nablau) = nablau
+ S(u, nablau) = nablau'
est = FeFunction(Vest, name = "est")
g = FeFunction(Vg, name = "g")
@@ -817,12 +849,8 @@ function denoise(ctx)
mesh = init_grid(g_arr, floor.(Int, size(g_arr) ./ 2^(n_refine / 2))...)
mesh_area = area(mesh)
- T(tdata, u) = u
- T(::typeof(adjoint), tdata, v) = v
- S(u, nablau) = u
-
st = L1L2TVState{1}(mesh;
- T, tdata = nothing, S,
+ T = Operators.TIdentity(), tdata = nothing, S = Operators.SIdentity(),
ctx.params.alpha1, ctx.params.alpha2,
ctx.params.lambda, ctx.params.beta,
ctx.params.gamma1, ctx.params.gamma2)
@@ -1276,33 +1304,37 @@ end
function inpaint(ctx)
# expect ctx.params.g_arr
+ #project_image! = project_l2_lagrange!
project_image! = project_qi_lagrange!
- n_refine = 5
# convert to cartesian coordinates
g_arr = from_img(ctx.params.g_arr)
mask_arr = from_img(ctx.params.mask_arr)
- mesh = init_grid(g_arr, floor.(Int, size(g_arr) ./ 2^(n_refine / 2))...)
+ mesh = init_grid(g_arr, floor.(Int, size(g_arr) ./ 2^(ctx.params.n_refine / 2))...)
mesh_area = area(mesh)
# inpaint specific stuff
Vg = FeSpace(mesh, P1(), (1,))
_mask = FeFunction(Vg, name = "mask")
- T(tdata, u) = abs(tdata[begin] - 1.) < 1e-8 ? u : zero(u)
- T(::typeof(adjoint), tdata, v) = T(tdata, v)
- S(u, nablau) = u
-
st = L1L2TVState{1}(mesh;
- T, tdata = _mask, S,
+ T = Operators.TInpaint(), tdata = _mask, ctx.params.S,
ctx.params.alpha1, ctx.params.alpha2,
ctx.params.lambda, ctx.params.beta,
ctx.params.gamma1, ctx.params.gamma2)
function interpolate_image_data!()
println("interpolate image data ...")
- project_image!(st.g, g_arr)
+ if ctx.params.n_refine == 0
+ # if we use a grid at image resolution, we wan't to avoid
+ # information loss at all cost.
+ # projet_l2_pixel is not suited for inpainting as it is a global
+ # operator, which may leak data from the inpainting domain
+ interpolate!(st.g, x -> evaluate_bilinear(g_arr, x))
+ else
+ project_image!(st.g, g_arr)
+ end
project_image!(st.tdata, mask_arr)
end
@@ -1313,6 +1345,7 @@ function inpaint(ctx)
pvd = paraview_collection(joinpath(ctx.outdir, "output.pvd")) do pvd
interpolate_image_data!()
+
pvd[0] = save_step(0)
i = 0
@@ -1338,7 +1371,7 @@ function inpaint(ctx)
# refinement stop criterion
k_refine += 1
- k_refine > n_refine && break
+ k_refine > ctx.params.n_refine && break
println("refine ...")
#estimate_res!(st)
marked_cells = Set(mark(st; theta = 0.5))
@@ -1362,11 +1395,11 @@ function inpaint(ctx)
#CSV.write(joinpath(ctx.outdir, "energies.csv"), df)
u_sampled = sample(st.u)
- saveimg(joinpath(ctx.outdir, "g.png"), to_img(g_arr))
+ saveimg(joinpath(ctx.outdir, "g.png"), grayclamp.(to_img(g_arr)))
saveimg(joinpath(ctx.outdir, "output.png"), grayclamp.(to_img(u_sampled)))
save_csv(joinpath(ctx.outdir, "mesh.csv"), st.u)
savedata(joinpath(ctx.outdir, "data.tex");
- ctx.params.eps_newton, n_refine,
+ ctx.params.eps_newton, ctx.params.n_refine,
st.alpha1, st.alpha2, st.lambda, st.beta, st.gamma1, st.gamma2,
width=size(u_sampled, 1), height=size(u_sampled, 2))
return st
@@ -1377,15 +1410,44 @@ function experiment_inpaint(ctx)
mask_arr = loadimg(joinpath(ctx.indir, "mask.png"))
mesh = init_grid(g_arr;)
- df = DataFrame()
- inpaint(Util.Context(ctx; name = "test", df,
+ inpaint(Util.Context(ctx; name = "test",
+ S = Operators.SIdentity,
+ n_refine = 5,
g_arr, mask_arr, mesh,
alpha1 = 0., alpha2 = 50., lambda = 1., beta = 1e-5,
gamma1 = 1e-4, gamma2 = 1e-4,
- eps_newton = 1e-4, adaptive = true,
+ eps_newton = 1e-4,
))
end
+function experiment_inpaint_denoise(ctx)
+ g_arr = loadimg(joinpath(ctx.indir, "input.png"))
+ mask_arr = loadimg(joinpath(ctx.indir, "mask.png"))
+ mesh = init_grid(g_arr;)
+
+ params = (
+ name = "test",
+ n_refine = 0,
+ g_arr, mask_arr, mesh,
+ noise_sigma = 0.1, noise_p = 0.02,
+ alpha1 = 0.2, alpha2 = 8., lambda = 1., #= beta see below =#
+ gamma1 = 1e-4, gamma2 = 1e-4,
+ eps_newton = 1e-4)
+
+ params.g_arr .= add_noise(g_arr; params.noise_sigma, params.noise_p)
+ saveimg(joinpath(ctx.outdir, "input_noisy.png"), grayclamp.(g_arr))
+
+ ctx(inpaint, "s-identity_beta-large"; params...,
+ S = Operators.SIdentity(), beta = 3e-1)
+ ctx(inpaint, "s-identity_beta-small"; params...,
+ S = Operators.SIdentity(), beta = 3e-4)
+
+ ctx(inpaint, "s-nabla_beta-large"; params...,
+ S = Operators.SNabla(), beta = 5e1)
+ ctx(inpaint, "s-nabla_beta-small"; params...,
+ S = Operators.SNabla(), beta = 5e-2)
+end
+
function optflow(ctx)
size(ctx.params.imgf0) == size(ctx.params.imgf1) ||
throw(ArgumentError("non-matching image sizes"))