22 Jun
2012
22 Jun
'12
6:40 p.m.
---
scribo/sandbox/z/skewness/skewness.cc | 128 +++++++++++++++++++++++++++++++++
1 files changed, 128 insertions(+), 0 deletions(-)
create mode 100644 scribo/sandbox/z/skewness/skewness.cc
diff --git a/scribo/sandbox/z/skewness/skewness.cc
b/scribo/sandbox/z/skewness/skewness.cc
new file mode 100644
index 0000000..0280fa3
--- /dev/null
+++ b/scribo/sandbox/z/skewness/skewness.cc
@@ -0,0 +1,128 @@
+#include <mln/io/pgm/all.hh>
+#include <mln/io/pbm/all.hh>
+#include <mln/core/image/dmorph/image_if.hh>
+#include <mln/core/image/dmorph/sub_image.hh>
+#include <mln/pw/all.hh>
+#include <mln/value/int_u8.hh>
+#include <mln/core/image/image2d.hh>
+#include <mln/core/neighb.hh>
+#include <mln/win/rectangle2d.hh>
+#include <mln/extension/adjust.hh>
+#include <mln/border/mirror.hh>
+#include <mln/extension/duplicate.hh>
+#include <mln/core/alias/window2d.hh>
+
+#include <mln/data/convert.hh>
+#include <mln/data/stretch.hh>
+#include <mln/debug/println.hh>
+#include <mln/data/transform.hh>
+#include <mln/fun/v2b/threshold.hh>
+#include <mln/io/dump/save.hh>
+
+mln::image2d<double> skewness;
+
+#include <scribo/binarization/sauvola.hh>
+
+namespace mln
+{
+
+ template <typename T, typename N>
+ std::pair<image2d<double>, image2d<double> >
+ init_(const image2d<T>& ima, const N& nbh)
+ {
+ image2d<double> output_2, output_3;
+ initialize(output_2, ima);
+ initialize(output_3, ima);
+ extension::adjust(ima, nbh);
+ border::mirror(ima);
+
+ typedef image2d<T> I;
+ mln_piter(I) p(ima.domain());
+ mln_niter(N) n(nbh, p);
+
+ extension::duplicate(ima);
+
+ for_all(p)
+ {
+ double sum = ima(p);
+ for_all(n)
+ sum += ima(n);
+
+ double mean = sum / (double)nbh.size();
+ double
+ res_2 = std::pow(ima(p) - mean, 2),
+ res_3 = std::pow(ima(p) - mean, 3);
+ for_all(n)
+ {
+ res_2 += std::pow(ima(n) - mean, 2);
+ res_3 += std::pow(ima(n) - mean, 3);
+ }
+
+ output_2(p) = res_2;
+ output_3(p) = res_3;
+ }
+
+ return std::pair<image2d<double>, image2d<double> >(output_2,
output_3);
+ }
+
+
+}
+
+
+
+int main(int argc, char *argv[])
+{
+ using namespace mln;
+
+ if (argc != 4)
+ {
+ std::cerr << "Usage: " << argv[0] << "
<input.pgm> <prefix> <win_size>" << std::endl;
+ return 1;
+ }
+
+ image2d<value::int_u8> input;
+ io::pgm::load(input, argv[1]);
+ std::string prefix = argv[2];
+
+ int win_size = atoi(argv[3]);
+ win::rectangle2d wr = win::rectangle2d(win_size,win_size);
+ window2d win;
+ win.insert(wr);
+
+ std::pair<image2d<double>, image2d<double> > res_init = init_(input,
neighb<window2d>(win));
+
+ image2d<double> x_mean_2 = res_init.first;
+ image2d<double> x_mean_3 = res_init.second;
+
+ {
+ io::pgm::save(data::stretch(value::int_u8(), x_mean_2), prefix +
"x_mean_2_stretch.pgm");
+ io::pgm::save(data::stretch(value::int_u8(), x_mean_3), prefix +
"x_mean_3_stretch.pgm");
+ }
+
+
+ initialize(skewness, x_mean_2);
+ data::fill(skewness, 0);
+
+ int n;
+ mln_piter_(image2d<double>) p(skewness.domain());
+ for_all(p)
+ {
+ n = win.size();
+ double frac = std::pow(std::sqrt((1.0/(double)n * x_mean_2(p))), 3);
+ skewness(p) = (std::sqrt(n * (n - 1)) / (n - 2)) * ((1.0/(double)n * x_mean_3(p)) /
frac);
+
+ // if (p == point2d(5, 5))
+ // std::cout << p << " - " << std::setprecision(8)
<< (double)x_mean_3(p) << " - " << (double)x_mean_2(p)
<< " - " << skewness(p) << std::endl;
+ }
+
+
+ {
+ io::dump::save(skewness, prefix + "skewness.dump");
+ io::pgm::save(data::convert(value::int_u8(), skewness), prefix +
"skewness_u8.pgm");
+ io::pbm::save(data::transform(skewness, fun::v2b::threshold<double>(0)), prefix
+ "skewness.pbm");
+ io::pgm::save(data::stretch(value::int_u8(), skewness), prefix +
"skewness_stretch.pgm");
+ }
+
+// image2d<bool> bin = scribo::binarization::sauvola(input, 51);
+// io::pbm::save(bin, "input.pbm");
+}
--
1.7.2.5