January 2013 - Olena-patches

branch exp/inpainting updated: olena-2.0-43-g39bb980
by git＠lrde.epita.fr 13 Jan '13

13 Jan '13

This is an automated email from the git hooks/post-receive script. It was generated because a ref change was pushed to the repository containing the project "Olena, a generic and efficient image processing platform". The branch exp/inpainting has been updated via 39bb98056d362c03aa92a0c83b667671ada0c3c6 (commit) from ca14b5cb20b525406712cff3249b4042ca57b573 (commit) Those revisions listed above that are new to this repository have not appeared on any other notification email; so we list those revisions in full, below. - Log ----------------------------------------------------------------- 39bb980 Add: TV L2 norm metric. ----------------------------------------------------------------------- Summary of changes: milena/mln/inpainting/derivative_kernels.hh | 41 +++++++++++++++++++ milena/mln/inpainting/metric/tv_l2_norm.hh | 57 +++++++++++++++++++++++++++ milena/mln/inpainting/tv.hh | 35 ++-------------- 3 files changed, 102 insertions(+), 31 deletions(-) create mode 100644 milena/mln/inpainting/metric/tv_l2_norm.hh hooks/post-receive -- Olena, a generic and efficient image processing platform

1 0

olena-2.0-43-g39bb980 Add: TV L2 norm metric.
by Coddy Levi 13 Jan '13

13 Jan '13

--- milena/mln/inpainting/derivative_kernels.hh | 41 +++++++++++++++++++ milena/mln/inpainting/metric/tv_l2_norm.hh | 57 +++++++++++++++++++++++++++ milena/mln/inpainting/tv.hh | 35 ++-------------- 3 files changed, 102 insertions(+), 31 deletions(-) create mode 100644 milena/mln/inpainting/metric/tv_l2_norm.hh diff --git a/milena/mln/inpainting/derivative_kernels.hh b/milena/mln/inpainting/derivative_kernels.hh index 04e1535..3e7bf63 100644 --- a/milena/mln/inpainting/derivative_kernels.hh +++ b/milena/mln/inpainting/derivative_kernels.hh @@ -1,6 +1,11 @@ #ifndef DERIVATIVE_KERNELS_HH_ # define DERIVATIVE_KERNELS_HH_ +#include <mln/core/alias/window2d.hh> +#include <mln/core/concept/iterator.hh> +#include <mln/convert/from_to.hh> +#include <mln/literal/zero.hh> + namespace mln { namespace inpainting @@ -97,7 +102,43 @@ namespace mln return ws; } + + namespace internal + { + template <typename I> + void convolve(const I& src, + I& ima, + const float* kernel) + { + static window2d full3x3; + static const bool vals[] = { 1, 1, 1, + 1, 1, 1, + 1, 1, 1 }; + + convert::from_to(vals, full3x3); + + mln_pixter(I) o(ima); + mln_pixter(const I) p(src); + mln_qixter(const I, const window2d) q(p, full3x3); + + for_all_2 (p, o) + { + mln_value(I) accu = literal::zero; + + unsigned k = 0; + for_all (q) + { + accu += q.val() * kernel[k]; + ++k; + } + + o.val() = accu; + } + } + + } } } + #endif diff --git a/milena/mln/inpainting/metric/tv_l2_norm.hh b/milena/mln/inpainting/metric/tv_l2_norm.hh new file mode 100644 index 0000000..17df267 --- /dev/null +++ b/milena/mln/inpainting/metric/tv_l2_norm.hh @@ -0,0 +1,57 @@ +#ifndef TV_L2_NORM_HH_ +# define TV_L2_NORM_HH_ + +# include <mln/core/concept/image.hh> +# include <mln/trait/value_.hh> +# include <mln/algebra/vec.hh> + +# include <mln/inpainting/derivative_kernels.hh> + +namespace mln +{ + namespace inpainting + { + namespace metric + { + template <typename I, typename M> + double tv_l2_norm(const Image<I>& ima_, + const Image<M>& mask_) + { + I dx; + I dy; + + const I& ima = exact(ima_); + + initialize(dx, ima); + initialize(dy, ima); + + internal::convolve(ima, dx, make_dx_kernel<float>()); + internal::convolve(ima, dy, make_dy_kernel<float>()); + + mln_piter(I) p( dx.domain() ); + + double sum = 0; + + const M& mask = exact(mask_); + + for_all(p) + { + if (mask(p)) + { + for (unsigned i = 0; i < mln_dim(mln_value(I)); ++i) + { + const float& a = dx(p)[i]; + const float& b = dy(p)[i]; + + sum += std::sqrt(a * a + b * b); + } + } + } + + return sum; + } + } + } +} + +#endif diff --git a/milena/mln/inpainting/tv.hh b/milena/mln/inpainting/tv.hh index 3b38270..8c8066b 100644 --- a/milena/mln/inpainting/tv.hh +++ b/milena/mln/inpainting/tv.hh @@ -7,6 +7,7 @@ #include <mln/convert/from_to.hh> #include <mln/norm/l1.hh> #include <mln/norm/l2.hh> +#include <mln/inpainting/metric/tv_l2_norm.hh> #include <mln/inpainting/derivative_kernels.hh> #include <mln/extension/adjust_duplicate.hh> @@ -78,37 +79,6 @@ namespace mln } } - template <typename I> - void convolve(const I& src, - I& ima, - const float* kernel) - { - static window2d full3x3; - static const bool vals[] = { 1, 1, 1, - 1, 1, 1, - 1, 1, 1 }; - - convert::from_to(vals, full3x3); - - mln_pixter(I) o(ima); - mln_pixter(const I) p(src); - mln_qixter(const I, const window2d) q(p, full3x3); - - for_all_2 (p, o) - { - mln_value(I) accu = literal::zero; - - unsigned k = 0; - for_all (q) - { - accu += q.val() * kernel[k]; - ++k; - } - - o.val() = accu; - } - } - // debug template <typename I, typename M> void print_first_value(const I& src, const M& mask) @@ -233,6 +203,9 @@ namespace mln data::paste(u, old); + std::cout << iterations << "\t" + << metric::tv_l2_norm(u, mask) << std::endl; + ++iterations; } -- 1.7.2.5

1 0

branch exp/inpainting updated: olena-2.0-42-gca14b5c
by git＠lrde.epita.fr 09 Jan '13

09 Jan '13

This is an automated email from the git hooks/post-receive script. It was generated because a ref change was pushed to the repository containing the project "Olena, a generic and efficient image processing platform". The branch exp/inpainting has been updated via ca14b5cb20b525406712cff3249b4042ca57b573 (commit) from 28a4ed0d52a219ebbc04a76d2bcb2b8f76ba25f9 (commit) Those revisions listed above that are new to this repository have not appeared on any other notification email; so we list those revisions in full, below. - Log ----------------------------------------------------------------- ca14b5c Fix output image encoding. Some cleaning. ----------------------------------------------------------------------- Summary of changes: milena/mln/inpainting/fast_isotropic_diffusion.hh | 3 -- milena/mln/inpainting/inpaint.hh | 22 ++++++++++++-------- milena/mln/inpainting/lpde.hh | 2 - milena/mln/inpainting/tv.hh | 10 +++----- 4 files changed, 17 insertions(+), 20 deletions(-) hooks/post-receive -- Olena, a generic and efficient image processing platform

1 0

olena-2.0-42-gca14b5c Fix output image encoding. Some cleaning.
by Coddy Levi 09 Jan '13

09 Jan '13

--- milena/mln/inpainting/fast_isotropic_diffusion.hh | 3 -- milena/mln/inpainting/inpaint.hh | 22 ++++++++++++-------- milena/mln/inpainting/lpde.hh | 2 - milena/mln/inpainting/tv.hh | 10 +++----- 4 files changed, 17 insertions(+), 20 deletions(-) diff --git a/milena/mln/inpainting/fast_isotropic_diffusion.hh b/milena/mln/inpainting/fast_isotropic_diffusion.hh index fad9fac..e28d787 100644 --- a/milena/mln/inpainting/fast_isotropic_diffusion.hh +++ b/milena/mln/inpainting/fast_isotropic_diffusion.hh @@ -79,9 +79,6 @@ namespace mln std::swap(u_, u_next_); } - - std::cout << iterations << "\t"; - } } /* mln::inpainting */ diff --git a/milena/mln/inpainting/inpaint.hh b/milena/mln/inpainting/inpaint.hh index ff58e75..fb99a96 100644 --- a/milena/mln/inpainting/inpaint.hh +++ b/milena/mln/inpainting/inpaint.hh @@ -127,9 +127,13 @@ namespace mln for_all_2(p, p_new) { - output(p).comp(0) = p_new.val()[0] * range; - output(p).comp(1) = p_new.val()[1] * range; - output(p).comp(2) = p_new.val()[2] * range; + algebra::vec<3u, int> v = p_new.val() * range; + + internal::interval(v, 0, 255); + + output(p).comp(0) = v[0]; + output(p).comp(1) = v[1]; + output(p).comp(2) = v[2]; } } @@ -142,12 +146,12 @@ namespace mln new_mask); if (out_error) - *out_error = color_error; - - std::cout << std::setprecision(6) << std::fixed - << t << "\t\t" - << (1 - color_error) * 100 << "\t\t" - << (1 - gradient_error) * 100 << "\t"; + *out_error = (1 - color_error) * 100; + else + std::cout << std::setprecision(6) << std::fixed + << t << "\t\t" + << (1 - color_error) * 100 << "\t\t" + << (1 - gradient_error) * 100 << "\t"; return output; } diff --git a/milena/mln/inpainting/lpde.hh b/milena/mln/inpainting/lpde.hh index 2eeca78..a293450 100644 --- a/milena/mln/inpainting/lpde.hh +++ b/milena/mln/inpainting/lpde.hh @@ -253,8 +253,6 @@ namespace mln ++t; } - - std::cout << t << "\t"; } void lpde::set(const std::vector<float_type>& w) diff --git a/milena/mln/inpainting/tv.hh b/milena/mln/inpainting/tv.hh index 44bf24b..3b38270 100644 --- a/milena/mln/inpainting/tv.hh +++ b/milena/mln/inpainting/tv.hh @@ -36,14 +36,14 @@ namespace mln namespace internal { - template <typename T> + template <typename T, typename V> inline void interval(T& v, - tv::float_type min, - tv::float_type max) + const V& min, + const V& max) { for (unsigned i = 0; i < mln_dim(T); ++i) - v[i] = std::max(std::min(v[i], max), min); + v[i] = std::max(std::min((V) v[i], max), min); } template <typename T> @@ -236,8 +236,6 @@ namespace mln ++iterations; } - std::cout << iterations << "\t"; - data::paste(u, src); } -- 1.7.2.5

1 0

olena-2.0-41-g28a4ed0 Add: out error variable to inpaint()
by Coddy Levi 07 Jan '13

07 Jan '13

--- milena/mln/inpainting/inpaint.hh | 9 +++++++-- 1 files changed, 7 insertions(+), 2 deletions(-) diff --git a/milena/mln/inpainting/inpaint.hh b/milena/mln/inpainting/inpaint.hh index d5e89c7..ff58e75 100644 --- a/milena/mln/inpainting/inpaint.hh +++ b/milena/mln/inpainting/inpaint.hh @@ -39,7 +39,8 @@ namespace mln I<T> inpaint(const Image<I<T> >& src_, const Image<I<bool> >& mask_, F inpaint_method, - unsigned inpaint_radius); + unsigned inpaint_radius, + double* out_error = NULL); /* * Implem @@ -49,7 +50,8 @@ namespace mln I<T> inpaint(const Image<I<T> >& inpaint_src_, const Image<I<bool> >& inpaint_mask_, F inpaint_method, - unsigned inpaint_radius) + unsigned inpaint_radius, + double* out_error) { const I<T>& inpaint_src = exact(inpaint_src_); const I<bool>& inpaint_mask = exact(inpaint_mask_); @@ -139,6 +141,9 @@ namespace mln new_src, new_mask); + if (out_error) + *out_error = color_error; + std::cout << std::setprecision(6) << std::fixed << t << "\t\t" << (1 - color_error) * 100 << "\t\t" -- 1.7.2.5

2 1

branch exp/inpainting updated: olena-2.0-41-g28a4ed0
by git＠lrde.epita.fr 07 Jan '13

07 Jan '13

This is an automated email from the git hooks/post-receive script. It was generated because a ref change was pushed to the repository containing the project "Olena, a generic and efficient image processing platform". The branch exp/inpainting has been updated via 28a4ed0d52a219ebbc04a76d2bcb2b8f76ba25f9 (commit) from 0885c7d8dd7b793359d5fc697aa9de65c90a08ff (commit) Those revisions listed above that are new to this repository have not appeared on any other notification email; so we list those revisions in full, below. - Log ----------------------------------------------------------------- 28a4ed0 Add: out error variable to inpaint() ----------------------------------------------------------------------- Summary of changes: milena/mln/inpainting/inpaint.hh | 9 +++++++-- 1 files changed, 7 insertions(+), 2 deletions(-) hooks/post-receive -- Olena, a generic and efficient image processing platform

1 0

branch exp/inpainting updated: olena-2.0-40-g0885c7d
by git＠lrde.epita.fr 07 Jan '13

07 Jan '13

This is an automated email from the git hooks/post-receive script. It was generated because a ref change was pushed to the repository containing the project "Olena, a generic and efficient image processing platform". The branch exp/inpainting has been updated via 0885c7d8dd7b793359d5fc697aa9de65c90a08ff (commit) via 8b1b95c4d94d77f584a6109f9c3a55e6dcc86200 (commit) via ee33739a3f9bf16c17ed0ca0e0b635be64e6b2cb (commit) via 818dc3790f71951840711e1031a7c965c2e45df8 (commit) via dd1d663316ca05dd1ea0af8f63eb6bb7a892a461 (commit) from 33248b9f4649130b71e97d3a15599a52502d9643 (commit) Those revisions listed above that are new to this repository have not appeared on any other notification email; so we list those revisions in full, below. - Log ----------------------------------------------------------------- 0885c7d Add: better integration of the tv diffusion term. cleaning. fix image copy. 8b1b95c Add: forward and backward difference kenrels for TV method, replaced tv method with primal dual version. ee33739 Add: initialize unknown area with mean value, value normalization, mean gradient error print. fix image copy. 818dc37 fix image copy dd1d663 Add: new mean gradient error metric. Some metric refactoring. ----------------------------------------------------------------------- Summary of changes: milena/mln/inpainting/derivative_kernels.hh | 91 ++++-- milena/mln/inpainting/fast_isotropic_diffusion.hh | 16 +- milena/mln/inpainting/inpaint.hh | 64 +++- milena/mln/inpainting/lpde.hh | 172 ++++----- .../mln/inpainting/metric/mean_gradient_error.hh | 75 ++++ milena/mln/inpainting/metric/mean_l1_error.hh | 57 +++ milena/mln/inpainting/metric/{psnr.hh => mse.hh} | 23 +- milena/mln/inpainting/metric/psnr.hh | 29 +-- milena/mln/inpainting/tv.hh | 372 +++++++++++++++++--- 9 files changed, 663 insertions(+), 236 deletions(-) create mode 100644 milena/mln/inpainting/metric/mean_gradient_error.hh create mode 100644 milena/mln/inpainting/metric/mean_l1_error.hh copy milena/mln/inpainting/metric/{psnr.hh => mse.hh} (60%) hooks/post-receive -- Olena, a generic and efficient image processing platform

1 0

olena-2.0-40-g0885c7d Add: better integration of the tv diffusion term. cleaning. fix image copy.
by Coddy Levi 07 Jan '13

07 Jan '13

--- milena/mln/inpainting/lpde.hh | 172 ++++++++++++++++++----------------------- 1 files changed, 76 insertions(+), 96 deletions(-) diff --git a/milena/mln/inpainting/lpde.hh b/milena/mln/inpainting/lpde.hh index f38b7a3..2eeca78 100644 --- a/milena/mln/inpainting/lpde.hh +++ b/milena/mln/inpainting/lpde.hh @@ -1,3 +1,6 @@ +#ifndef LPDE_HH_ +# define LPDE_HH_ + #include <algorithm> #include <vector> #include <limits> @@ -67,6 +70,8 @@ namespace mln this->dt = dt; this->dx = dx; this->t_f = t_f; + + this->w = std::vector<float_type>(6, 0.f); } template <typename T> @@ -114,9 +119,7 @@ namespace mln I& u_next, const std::vector<lpde::float_type>& w, const M& mask, - lpde::float_type& max_diff_norm, - I& u_x, - I& u_y) + lpde::float_type& max_diff_norm) { typedef lpde::float_type float_type; typedef mln_value(I) T; @@ -143,86 +146,78 @@ namespace mln mln_pixter(I) p_next(u_next); mln_pixter(const M) p_m(mask); - mln_pixter(I) pu_x(u_x); - mln_pixter(I) pu_y(u_y); - mln_qixter(const I, const window2d) q(p, full3x3); float_type sse = 0; - for (p.start(), p_next.start(), p_m.start(), - pu_x.start(), pu_y.start(); + for (p.start(), p_next.start(), p_m.start(); p.is_valid(); - p.next(), p_next.next(), p_m.next(), - pu_x.next(), pu_y.next()) + p.next(), p_next.next(), p_m.next()) { if (p_m.val()) - { - const T old = p.val(); - - if (fabs(w[1]) > activation_threshold) - p_next.val() = p.val() * w[1]; - - T p_derivative[derivative_kernels_size]; - - for (unsigned i = 0; i < derivative_kernels_size; ++i) - p_derivative[i] = literal::zero; - - unsigned k = 0; - for_all (q) { - for (unsigned i = 0; i < derivative_kernels_size; ++i) - { - p_derivative[i] += q.val() * derivative_kernel[i][k]; - } - ++k; - } - - const T& p_x = p_derivative[0]; - const T& p_y = p_derivative[1]; - const T& p_xx = p_derivative[2]; - const T& p_yy = p_derivative[3]; - const T& p_xy = p_derivative[4]; - - const T p_x2 = components_product(p_x, p_x); - const T p_y2 = components_product(p_y, p_y); - const T p_xx2 = components_product(p_xx, p_xx); - const T p_yy2 = components_product(p_yy, p_yy); - const T p_xy2 = components_product(p_xy, p_xy); + const T old = p.val(); - pu_x.val() = normalize(p_x); - pu_y.val() = normalize(p_y); + if (fabs(w[1]) > activation_threshold) + p_next.val() = p.val() * w[1]; - //||∇u||^2 - if (fabs(w[2]) > activation_threshold) - p_next.val() += w[2] * (p_x2 + p_y2); + T p_derivative[derivative_kernels_size]; - // laplacian - if (fabs(w[3]) > activation_threshold) - p_next.val() += w[3] * (p_xx + p_yy); - - // tr(H^2(u)) - if (fabs(w[4]) > activation_threshold) - p_next.val() += w[4] * (p_xx2 + 2 * p_xy2 + p_yy2); - - // (∇u)^T H_u (∇u) - if (fabs(w[5]) > activation_threshold) - p_next.val() += w[5] * (components_product(p_x2, p_xx) - + 2 * components_product(components_product(p_x, p_y), p_xy) - + components_product(p_y2, p_yy)); - - float_type diff = norm::l2(old - p_next.val()); - sse += diff; + for (unsigned i = 0; i < derivative_kernels_size; ++i) + p_derivative[i] = literal::zero; - if (diff > max_diff_norm) + unsigned k = 0; + for_all (q) { - max_diff_norm = diff; + for (unsigned i = 0; i < derivative_kernels_size; ++i) + { + p_derivative[i] += q.val() * derivative_kernel[i][k]; + } + ++k; } - } + + const T& p_x = p_derivative[0]; + const T& p_y = p_derivative[1]; + const T& p_xx = p_derivative[2]; + const T& p_yy = p_derivative[3]; + const T& p_xy = p_derivative[4]; + + const T p_x2 = components_product(p_x, p_x); + const T p_y2 = components_product(p_y, p_y); + const T p_xx2 = components_product(p_xx, p_xx); + const T p_yy2 = components_product(p_yy, p_yy); + const T p_xy2 = components_product(p_xy, p_xy); + + // //||∇u||^2 + if (fabs(w[2]) > activation_threshold) + p_next.val() += w[2] * (p_x2 + p_y2); + + // // laplacian + if (fabs(w[3]) > activation_threshold) + p_next.val() += w[3] * (p_xx + p_yy); + + // // tr(H^2(u)) + if (fabs(w[4]) > activation_threshold) + p_next.val() += w[4] * (p_xx2 + 2 * p_xy2 + p_yy2); + + // // (∇u)^T H_u (∇u) + if (fabs(w[5]) > activation_threshold) + p_next.val() += w[5] * (components_product(p_x2, p_xx) + + 2 * components_product(components_product(p_x, p_y), p_xy) + + components_product(p_y2, p_yy)); + + internal::interval(p_next.val(), 0.0, 1.0); + + float_type diff = norm::l2(old - p_next.val()); + sse += diff; + + if (diff > max_diff_norm) + { + max_diff_norm = diff; + } + } } - - inpainting::tv(u_next, u_x, u_y, mask, max_diff_norm); - + return sse; } @@ -230,14 +225,8 @@ namespace mln void lpde::operator()(I<T>& src, const I<bool>& mask) { - // mean_diffusion mean_diff; - // mean_diff(src, mask); - - // fast_isotropic_diffusion fast_iso; - // fast_iso(src, mask); - I<T>& u_ = src; - I<T> u_next_(src); + I<T> u_next_ = duplicate(src); extension::adjust_duplicate(u_, 1); extension::adjust_duplicate(u_next_, 1); @@ -245,36 +234,24 @@ namespace mln I<T>* u = &u_; I<T>* u_next = &u_next_; - float_type sse = std::numeric_limits<float_type>::max() / 2.0; - - float_type sse_old = std::numeric_limits<float_type>::max(); - float_type max_diff_norm = std::numeric_limits<float_type>::max(); - float_type min_max_diff_norm = std::numeric_limits<float_type>::max(); + inpainting::tv_step<I, T> tv(src, mask); - I<T> u_x(u_.domain()); - I<T> u_y(u_.domain()); - - float t = 0; + unsigned t = 0; while (t < this->t_f && max_diff_norm > CONVERGENCE_THRESHOLD - && !std::isnan(max_diff_norm) - ) + && !std::isnan(max_diff_norm)) { - max_diff_norm = 0; - sse_old = sse; - - sse = diffuse(*u, *u_next, w, mask, max_diff_norm, u_x, u_y); - - if (max_diff_norm < min_max_diff_norm) - min_max_diff_norm = max_diff_norm; - + float sse = diffuse(*u, *u_next, w, mask, max_diff_norm); + + tv(*u_next); + std::swap(u, u_next); - t += this->dt; + ++t; } std::cout << t << "\t"; @@ -282,8 +259,11 @@ namespace mln void lpde::set(const std::vector<float_type>& w) { - this->w = w; + for (unsigned i = 0; i < std::min(w.size(), this->w.size()); ++i) + this->w[i] = w[i]; } } /* mln::inpainting */ } /* mln */ + +#endif -- 1.7.2.5

1 0

olena-2.0-39-g8b1b95c Add: forward and backward difference kenrels for TV method, replaced tv method with primal dual version.
by Coddy Levi 07 Jan '13

07 Jan '13

--- milena/mln/inpainting/derivative_kernels.hh | 91 +++++-- milena/mln/inpainting/tv.hh | 372 ++++++++++++++++++++++---- 2 files changed, 381 insertions(+), 82 deletions(-) diff --git a/milena/mln/inpainting/derivative_kernels.hh b/milena/mln/inpainting/derivative_kernels.hh index a7f3c34..04e1535 100644 --- a/milena/mln/inpainting/derivative_kernels.hh +++ b/milena/mln/inpainting/derivative_kernels.hh @@ -1,30 +1,71 @@ +#ifndef DERIVATIVE_KERNELS_HH_ +# define DERIVATIVE_KERNELS_HH_ + namespace mln { namespace inpainting { - namespace - { - const float a = 61; - const float b = 17; - - template <typename T> - const T* make_dx_kernel() - { - static const T ws[] = { b, a, b, - 0, 0, 0, - -b, -a, -b }; - return ws; - } - - template <typename T> - const T* make_dy_kernel() - { - static const T ws[] = { b, 0, -b, - a, 0, -a, - b, 0, -b }; + namespace internal + { + const float a = 2; + const float b = 1; + } + + template <typename T> + const T* make_dx_kernel() + { + static const T ws[] = { internal::b, internal::a, internal::b, + 0, 0, 0, + -internal::b, -internal::a, -internal::b }; + return ws; + } + + template <typename T> + const T* make_dy_kernel() + { + static const T ws[] = { internal::b, 0, -internal::b, + internal::a, 0, -internal::a, + internal::b, 0, -internal::b }; - return ws; - } + return ws; + } + + template <typename T> + const T* make_forward_horizontal_kernel() + { + static const T ws[] = { 0, 0, 0, + 0, -1, 1, + 0, 0, 0 }; + return ws; + } + + template <typename T> + const T* make_forward_vertical_kernel() + { + static const T ws[] = { 0, 0, 0, + 0, -1, 0, + 0, 1, 0 }; + + return ws; + } + + template <typename T> + const T* make_backward_horizontal_kernel() + { + static const T ws[] = { 0, 0, 0, + -1, 1, 0, + 0, 0, 0 }; + return ws; + } + + template <typename T> + const T* make_backward_vertical_kernel() + { + static const T ws[] = { 0, -1, 0, + 0, 1, 0, + 0, 0, 0 }; + + return ws; } template <typename T> @@ -51,10 +92,12 @@ namespace mln const T* make_dxy_kernel() { static const T ws[] = {0.1051, 0, -0.1051, - 0, 0, 0, - -0.1051, 0, 0.1051}; + 0, 0, 0, + -0.1051, 0, 0.1051}; return ws; } } } + +#endif diff --git a/milena/mln/inpainting/tv.hh b/milena/mln/inpainting/tv.hh index d5bed23..44bf24b 100644 --- a/milena/mln/inpainting/tv.hh +++ b/milena/mln/inpainting/tv.hh @@ -1,8 +1,16 @@ +#ifndef TV_HH_ +# define TV_HH_ + #include <mln/core/alias/window2d.hh> #include <mln/core/concept/iterator.hh> #include <mln/algebra/vec.hh> #include <mln/convert/from_to.hh> +#include <mln/norm/l1.hh> #include <mln/norm/l2.hh> +#include <mln/inpainting/derivative_kernels.hh> + +#include <mln/extension/adjust_duplicate.hh> +#include <mln/inpainting/fast_isotropic_diffusion.hh> namespace mln { @@ -10,84 +18,332 @@ namespace mln { struct tv { + public: + tv(unsigned max_iterations = 10000) + : max_iterations_ (max_iterations) + { + } + typedef float float_type; + + template <template <typename> class I, typename T> + void operator()(I<T>& src, + const I<bool>& mask); + + private: + unsigned max_iterations_; }; - template <typename T> - inline - void interval(T& v, - tv::float_type min, - tv::float_type max) + namespace internal { - for (unsigned i = 0; i < mln_dim(T); ++i) - v[i] = std::max(std::min(v[i], max), min); + template <typename T> + inline + void interval(T& v, + tv::float_type min, + tv::float_type max) + { + for (unsigned i = 0; i < mln_dim(T); ++i) + v[i] = std::max(std::min(v[i], max), min); + } + + template <typename T> + inline + void normalize_linf(T& u, + T& v) + { + for (unsigned i = 0; i < mln_dim(T); ++i) + { + float norm = std::sqrt(u[i] * u[i] + v[i] * v[i]); + + norm = std::max(1.0f, norm); + + u[i] /= norm; + v[i] /= norm; + } + } + + template <typename I, typename T> + void linear(I& y, + const I& x, + const T& a, + const T& b, + const mln_value(I)& c = literal::zero) + { + mln_pixter(const I) p_x(x); + mln_pixter(I) p_y(y); + + for_all_2 (p_x, p_y) + { + p_y.val() = a * p_y.val() + b * p_x.val() + c; + } + } + + template <typename I> + void convolve(const I& src, + I& ima, + const float* kernel) + { + static window2d full3x3; + static const bool vals[] = { 1, 1, 1, + 1, 1, 1, + 1, 1, 1 }; + + convert::from_to(vals, full3x3); + + mln_pixter(I) o(ima); + mln_pixter(const I) p(src); + mln_qixter(const I, const window2d) q(p, full3x3); + + for_all_2 (p, o) + { + mln_value(I) accu = literal::zero; + + unsigned k = 0; + for_all (q) + { + accu += q.val() * kernel[k]; + ++k; + } + + o.val() = accu; + } + } + + // debug + template <typename I, typename M> + void print_first_value(const I& src, const M& mask) + { + for (unsigned i = 0; i < src.nelements(); ++i) + { + point2d p = src.point_at_index(i); + + if (mask(p)) + { + std::cout << p << std::endl; + std::cout << src(p) << std::endl; + return; + } + } + } + + template <typename I> + float max_diff_norm(const I& u, const I& old) + { + float max = 0; + + mln_pixter(const I) p(u); + mln_pixter(const I) p_old(old); + + for_all_2 (p, p_old) + { + max = std::max(max, norm::l1(p.val() - p_old.val())); + } + + return max * max; // l2 norm + } } - template <typename I, typename M> - void tv(I& u, - const I& u_x, - const I& u_y, - const M& mask, - tv::float_type& max_diff_norm) + // TV regularization inpainting with PrimalDual algorithm + template <template <typename> class I, typename T> + void tv::operator()(I<T>& src, + const I<bool>& mask) { - typedef tv::float_type float_type; + typedef float float_type; - extension::adjust_duplicate(u_x, 1); - extension::adjust_duplicate(u_y, 1); + const float_type* dx_bkern = make_backward_horizontal_kernel<float_type>(); + const float_type* dy_bkern = make_backward_vertical_kernel<float_type>(); + const float_type* dx_fkern = make_forward_horizontal_kernel<float_type>(); + const float_type* dy_fkern = make_forward_vertical_kernel<float_type>(); - const float_type* dx_kern = make_dx_kernel<float_type>(); - const float_type* dy_kern = make_dy_kernel<float_type>(); + const float_type sigma = 1.0f / std::sqrt(8.0f); + + I<T> u_x; + initialize(u_x, src); - static window2d full3x3; - static const bool vals[] = { 1, 1, 1, - 1, 1, 1, - 1, 1, 1 }; + I<T> u_y; + initialize(u_y, src); - convert::from_to(vals, full3x3); + I<T> u_x2; + initialize(u_x2, src); + I<T> u_y2; + initialize(u_y2, src); + + I<T> u = duplicate(src); + + I<T> old = duplicate(u); + + I<T> u_aux = duplicate(u); + + internal::convolve(u, u_x, dx_fkern); + internal::convolve(u, u_y, dy_fkern); + + float_type max_diff = std::numeric_limits<float_type>::max(); + + unsigned iterations = 0; + + while ((max_diff > CONVERGENCE_THRESHOLD) + && (iterations < this->max_iterations_)) + { + internal::convolve(u_aux, u_x2, dx_fkern); + internal::convolve(u_aux, u_y2, dy_fkern); + + internal::linear(u_x, u_x2, 1.0f, sigma); + internal::linear(u_y, u_y2, 1.0f, sigma); + + mln_pixter(I<T>) p_x(u_x); + mln_pixter(I<T>) p_y(u_y); + + for_all_2 (p_x, p_y) + { + internal::normalize_linf(p_x.val(), p_y.val()); + } + + // divergence + { + internal::convolve(u_x, u_x2, dx_bkern); + internal::convolve(u_y, u_y2, dy_bkern); + + internal::linear(u, u_x2, 0.f, 1.f); + internal::linear(u, u_y2, 1.f, 1.f); + } - mln_pixter(I) p(u); - mln_pixter(const M) p_m(mask); - mln_pixter(const I) p_x(u_x); - mln_pixter(const I) p_y(u_y); + internal::linear(u, u, sigma, 0.0f); + internal::linear(u, old, 1.0f, 1.0f); - mln_qixter(const I, const window2d) q_x(p_x, full3x3); - mln_qixter(const I, const window2d) q_y(p_y, full3x3); + mln_pixter(I<T>) p(u); + mln_pixter(const I<T>) p_src(src); + mln_pixter(const I<bool>) p_m(mask); - for (p.start(), p_m.start(), p_x.start(), p_y.start(); - p.is_valid(); - p.next(), p_m.next(), p_x.next(), p_y.next()) + for_all_3 (p, p_src, p_m) + { + if (!p_m.val()) + { + p.val() = p_src.val(); + } + else + { + internal::interval(p.val(), 0.0f, 1.0f); + } + } + + internal::linear(u_aux, u, 0.0f, 2.0f); + internal::linear(u_aux, old, 1.0f, -1.0f); + + max_diff = internal::max_diff_norm(u, old); + + data::paste(u, old); + + ++iterations; + } + + std::cout << iterations << "\t"; + + data::paste(u, src); + } + + template <template <typename> class I, typename T> + struct tv_step + { + public: + typedef float float_type; + + tv_step(I<T>& src, + const I<bool>& mask_) + : mask (mask_) + { + this->src = duplicate(src); + + initialize(u_x, src); + initialize(u_y, src); + initialize(u_x2, src); + initialize(u_y2, src); + + u = duplicate(src); + u_aux = duplicate(u); + + dx_bkern = make_backward_horizontal_kernel<float_type>(); + dy_bkern = make_backward_vertical_kernel<float_type>(); + dx_fkern = make_forward_horizontal_kernel<float_type>(); + dy_fkern = make_forward_vertical_kernel<float_type>(); + + sigma = 1.0f / std::sqrt(8.0f); + } + + void operator()(I<T>& ima) + { + internal::convolve(u_aux, u_x2, dx_fkern); + internal::convolve(u_aux, u_y2, dy_fkern); + + internal::linear(u_x, u_x2, 1.0f, sigma); + internal::linear(u_y, u_y2, 1.0f, sigma); + + mln_pixter(I<T>) p_x(u_x); + mln_pixter(I<T>) p_y(u_y); + + for_all_2 (p_x, p_y) + { + internal::normalize_linf(p_x.val(), p_y.val()); + } + + // divergence { - if (p_m.val()) - { - mln_value(I) old = p.val(); + internal::convolve(u_x, u_x2, dx_bkern); + internal::convolve(u_y, u_y2, dy_bkern); + + internal::linear(u, u_x2, 0.f, 1.f); + internal::linear(u, u_y2, 1.f, 1.f); + } - mln_value(I) p_x = literal::zero; - mln_value(I) p_y = literal::zero; + internal::linear(u, u, sigma, 0.0f); + internal::linear(u, ima, 1.0f, 1.0f); - unsigned k = 0; - for_all_2 (q_x, q_y) + mln_pixter(I<T>) p(u); + mln_pixter(const I<T>) p_src(src); + mln_pixter(const I<bool>) p_m(mask); + + for_all_3 (p, p_src, p_m) + { + if (!p_m.val()) { - p_x += q_x.val() * dx_kern[k]; - p_y += q_y.val() * dy_kern[k]; - ++k; + p.val() = p_src.val(); } - - mln_value(I) dv = p_x + p_y; - - p.val() += dv; - - interval(p.val(), 0.0, 255.0); - - float_type diff = norm::l2(old - p.val()); - - if (diff > max_diff_norm) - { - max_diff_norm = diff; - } - } - } - } + else + { + internal::interval(p.val(), 0.0f, 1.0f); + } + } + + internal::linear(u_aux, u, 0.0f, 2.0f); + internal::linear(u_aux, ima, 1.0f, -1.0f); + + float max_diff = internal::max_diff_norm(u, ima); + + data::paste(u, ima); + } + + private: + I<T> src; + const I<bool>& mask; + + I<T> u_x; + I<T> u_y; + + I<T> u_x2; + I<T> u_y2; + + I<T> u; + I<T> u_aux; + + const float_type* dx_bkern; + const float_type* dy_bkern; + const float_type* dx_fkern; + const float_type* dy_fkern; + + float_type sigma; + }; } /* mln::inpainting */ } /* mln */ + +#endif -- 1.7.2.5

1 0

olena-2.0-38-gee33739 Add: initialize unknown area with mean value, value normalization, mean gradient error print. fix image copy.
by Coddy Levi 07 Jan '13

07 Jan '13

--- milena/mln/inpainting/inpaint.hh | 64 +++++++++++++++++++++++++++++-------- 1 files changed, 50 insertions(+), 14 deletions(-) diff --git a/milena/mln/inpainting/inpaint.hh b/milena/mln/inpainting/inpaint.hh index 9db2438..d5e89c7 100644 --- a/milena/mln/inpainting/inpaint.hh +++ b/milena/mln/inpainting/inpaint.hh @@ -1,3 +1,6 @@ +#ifndef INPAINT_HH_ +# define INPAINT_HH_ + #include <mln/core/image/dmorph/image_if.hh> #include <mln/data/fill.hh> @@ -17,6 +20,13 @@ #include <mln/inpainting/metric/psnr.hh> +#include <mln/border/resize.hh> + +#include <mln/inpainting/metric/mean_gradient_error.hh> +#include <mln/inpainting/metric/mean_l1_error.hh> + +#include <mln/extension/adjust_duplicate.hh> + namespace mln { namespace inpainting @@ -63,9 +73,10 @@ namespace mln roi.crop_wrt(inpaint_src.domain()); } - //typedef mln_sum(T) sum_type; typedef algebra::vec<3u, typename F::float_type> sum_type; + float range = 255.0; + I<sum_type> new_src(roi); { // paste inpaint_src to new_src mln_piter(I<T>) p(roi); @@ -73,16 +84,31 @@ namespace mln for_all_2(p, p_new) { - p_new.val()[0] = inpaint_src(p).comp(0); - p_new.val()[1] = inpaint_src(p).comp(1); - p_new.val()[2] = inpaint_src(p).comp(2); + p_new.val()[0] = inpaint_src(p).comp(0) / range; + p_new.val()[1] = inpaint_src(p).comp(1) / range; + p_new.val()[2] = inpaint_src(p).comp(2) / range; } } + I<sum_type> old_src = duplicate(new_src); + I<bool> new_mask(roi); data::paste(inpaint_mask | roi, new_mask); - data::fill((new_src | pw::value(new_mask)).rw(), literal::zero); + sum_type v; + { + typedef typename trait::value_<T>::comp_0 inner_type; + inner_type middle = mln_max(inner_type) / 2; + + for (unsigned i = 0; i < mln_dim(sum_type); ++i) + v[i] = middle; + } + + data::fill((new_src | pw::value(new_mask)).rw(), v / range); + + border::resize(new_src, 1); + extension::adjust_duplicate(new_src, 1); + border::resize(new_mask, 0); mln::util::timer t; t.start(); @@ -91,24 +117,32 @@ namespace mln t.stop(); - double psnr = metric::psnr(inpaint_src | pw::value(inpaint_mask), - new_src | pw::value(new_mask)); + I<T> output = duplicate(inpaint_src); - std::cout << t << "\t" - << psnr << "\t"; - - I<T> output(inpaint_src); { // paste new_src to output mln_piter(I<T>) p(roi); mln_pixter(const I<sum_type>) p_new(new_src); for_all_2(p, p_new) { - output(p).comp(0) = p_new.val()[0]; - output(p).comp(1) = p_new.val()[1]; - output(p).comp(2) = p_new.val()[2]; + output(p).comp(0) = p_new.val()[0] * range; + output(p).comp(1) = p_new.val()[1] * range; + output(p).comp(2) = p_new.val()[2] * range; } } + + float gradient_error = metric::mean_gradient_error(old_src, + new_src, + new_mask); + + float color_error = metric::mean_l1_error(old_src, + new_src, + new_mask); + + std::cout << std::setprecision(6) << std::fixed + << t << "\t\t" + << (1 - color_error) * 100 << "\t\t" + << (1 - gradient_error) * 100 << "\t"; return output; } @@ -116,3 +150,5 @@ namespace mln } /* mln::inpainting */ } /* mln*/ + +#endif -- 1.7.2.5

1 0