add test...

1 files changed, 618 insertions, 0 deletions

diff --git a/meowpp.test/src/autostitch.cpp b/meowpp.test/src/autostitch.cpp
new file mode 100644
index 0000000..e63ef43
--- /dev/null
+++ b/meowpp.test/src/autostitch.cpp
@@ -0,0 +1,618 @@
+#include <cstdio>
+
+#include "autostitch.h"
+
+#include <opencv/cv.h>
+#include <opencv/highgui.h>
+
+#include "meowpp/Usage.h"
+
+#include "meowpp/colors/RGB_Space.h"
+
+#include "meowpp/dsa/DisjointSet.h"
+
+#include "meowpp/geo/Vectors.h"
+
+#include "meowpp/gra/Bitmap.h"
+#include "meowpp/gra/Photo.h"
+#include "meowpp/gra/Camera.h"
+#include "meowpp/gra/WatchBall.h"
+
+#include "meowpp/math/utility.h"
+#include "meowpp/math/methods.h"
+
+
+extern "C"{
+#include <sys/types.h>
+#include <dirent.h>
+}
+
+#include <vector>
+#include <algorithm>
+#include <string>
+#include <cstdlib>
+
+
+using namespace meow;
+
+//////////////////////////////////////////////////////////////////////
+
+Usage usg("autostitch");
+
+double p0 = 0.07, P = 0.99;
+double q = 0.7, r = 0.01, Q = 0.97;
+double stop = 1;
+double o_radius = 500;
+
+MyK_Match match;
+std::vector<Bitmap<RGBf_Space> >  input_bitmap;
+std::vector<Bitmap<RGBf_Space> > output_bitmap;
+
+std::vector<std::vector<FeaturePoint<double, double> > > fps;
+std::vector<std::vector<Vector<double>               > > fpsv;
+std::vector<std::vector<FeaturePointIndexPairs       > > pairs;
+
+struct OutputSet {
+  struct Edge {
+    std::vector<Vector<double> > v1;
+    std::vector<Vector<double> > v2;
+    size_t i1, i2;
+    bool done;
+    Edge(size_t ii1, size_t ii2): i1(ii1), i2(ii2) {
+      done = false;
+    }
+    bool operator<(Edge const& e) const {
+      return (v1.size() < e.v1.size());
+    }
+  };
+  std::vector<Camera<RGBf_Space> > cameras;
+  std::vector<Edge>                edges;
+};
+
+std::vector<OutputSet> outputs;
+
+//////////////////////////// **# setup #** ///////////////////////////
+bool setup(int argc, char** argv) {
+  usg.optionAdd('h', "Display this help document.");
+  usg.optionAdd('i',
+                "Specify the input images are in <type> "
+                "instead of specifying from arguments",
+                "<dirname>", "",
+                false);
+  usg.optionAdd('o',
+                "Output file name, (not include '.jpg' suffix)",
+                "<filename>",
+                "output",
+                false);
+  usg.optionAdd('d',
+                "Specify which Feature-Point-Detect algorithm to use",
+                "<algorithm>",
+                "",
+                true);
+  usg.optionAdd('p',
+                "Pribabilicity for RANSAC to choose a right feature point",
+                "<floating point>", stringPrintf("%.10f", p0),
+                false);
+  usg.optionAdd('P',
+                "Pribabilicity for RANSAC access",
+                "<floating point>", stringPrintf("%.10f", P),
+                false);
+  usg.optionAdd('q',
+                "p1 for Prob. Model",
+                "<floationg Point>", stringPrintf("%.10f", q),
+                false);
+  usg.optionAdd('r',
+                "p0 for Prob. Model",
+                "<floationg Point>", stringPrintf("%.10f", r),
+                false);
+  usg.optionAdd('Q',
+                "p_min for Prob. Model",
+                "<floationg Point>", stringPrintf("%.10f", Q),
+                false);
+  usg.optionAdd('s',
+                "stop threshold for boundle adjustment",
+                "<floationg Point>", stringPrintf("%.10f", stop),
+                false);
+  usg.optionAdd('O',
+                "output ball radius",
+                "<floationg Point>", stringPrintf("%.10f", o_radius),
+                false);
+  std::vector<std::string> fpsd_algorithm_list = ObjSelector<FPSD_ID>::names();
+  for (size_t i = 0, I = fpsd_algorithm_list.size(); i < I; i++) {
+    const ObjBase* tmp = ObjSelector<FPSD_ID>::get(fpsd_algorithm_list[i]);
+    usg.optionValueAcceptAdd('d',
+                             fpsd_algorithm_list[i],
+                             tmp->type());
+    usg.import(((MyFeaturePointsDetector*)tmp)->usage());
+  }
+  usg.import(match.usage());
+  usg.import(MyRansacCheck::usage());
+  // set arg
+  std::string err_msg;
+  bool ok = usg.arguments(argc, argv, &err_msg);
+  if (usg.hasOptionSetup('h')) {
+    printf("%s\n", usg.usage().c_str());
+    exit(0);
+  }
+  if (!ok) {
+    fprintf(stderr, "%s\n", err_msg.c_str());
+    exit(-1);
+  }
+  return true;
+}
+
+
+//////////////// **# Input images and convert it #** /////////////////
+bool input() {
+  std::vector<std::string> input_name;
+  if (!usg.hasOptionSetup('i')) {
+    input_name = usg.procArgs();
+  }
+  else {
+    std::string base = usg.optionValue('i', 0);
+    if (base.length() == 0 || base[base.length() - 1] != '/') {
+      base += "/";
+    }
+    DIR* dir = opendir(base.c_str());
+    if (!dir) {
+      fprintf(stderr, "can't open dir '%s'\n", base.c_str());
+      return -1;
+    }
+    for (dirent* ent; (ent = readdir(dir)) != NULL; ) {
+      if (!cstringEndWith(ent->d_name, 4, ".jpeg", ".jpg", ".JPG", ".JPEG")) {
+        continue;
+      }
+      input_name.push_back(base + std::string(ent->d_name));
+    }
+  }
+  messagePrintf(1, "Loading images");
+  for (size_t i = 0; i < input_name.size(); i++) {
+    messagePrintf(1, "%s", input_name[i].c_str());
+    cv::Mat img = cv::imread(input_name[i], CV_LOAD_IMAGE_COLOR);
+    if (!img.data) {
+      messagePrintf(-1, "opencv read error!, ignore");
+      continue;
+    }
+    size_t width  = img.size().width ;
+    size_t height = img.size().height;
+    size_t index = input_bitmap.size();
+    input_bitmap.resize(index + 1);
+    input_bitmap[index].size(height, width, RGBf_Space(0));
+    for (size_t x = 0; x < width; x++) {
+      for (size_t y = 0; y < height; y++) {
+        RGBi_Space tmp(Vector3D<int>(
+            img.at<cv::Vec3b>(y, x)[2],
+            img.at<cv::Vec3b>(y, x)[1],
+            img.at<cv::Vec3b>(y, x)[0]));
+        RGBf_Space p;
+        colorTransformate(tmp, &p);
+        input_bitmap[index].pixel(y, x, p);
+      }
+    }
+    messagePrintf(-1, "%lux%lu, ok", width, height);
+  }
+  messagePrintf(-1, "ok");
+  return true;
+}
+
+//////////////////////// **# FeaturePoint #** ////////////////////////
+bool detect() {
+  std::string fpsd_algo_name = usg.optionValue('d', 0);
+  MyFeaturePointsDetector* detector(
+    (MyFeaturePointsDetector*)ObjSelector<FPSD_ID>::create(fpsd_algo_name));
+  detector->usage(usg);
+  fps .resize(input_bitmap.size());
+  fpsv.resize(input_bitmap.size());
+  for (size_t i = 0, I = input_bitmap.size(); i < I; i++) {
+    messagePrintf(1, "Detect the feature points for %lu-th pic", i);
+    fps[i] = detector->detect(input_bitmap[i]);
+    messagePrintf(-1, "ok, %lu", fps[i].size());
+    for (size_t j = 0, J = fps[i].size(); j < J; j++) {
+      fpsv[i].push_back(fps[i][j].position());
+    }
+  }
+  delete detector;
+  return true;
+}
+
+
+//////////////////////////// **# k-match #** /////////////////////////
+bool kmatch() {
+  match.usage(usg);
+  messagePrintf( 1, "run k-match");
+  FeaturePointIndexPairs mat(match.match(fps));
+  pairs.resize(input_bitmap.size());
+  for (size_t i = 0, I = input_bitmap.size(); i < I; i++) {
+    pairs[i].resize(I);
+  }
+  for (size_t i = 0, I = mat.size(); i < I; ++i) {
+    pairs[mat[i].from.first][mat[i].to.first].push_back(mat[i]);
+  }
+  messagePrintf(-1, "ok");
+  return true;
+}
+
+//////////////////////////// **# RANSAC #** //////////////////////////
+bool ransac() {
+  messagePrintf( 1, "RANSAC");
+  MyRansacCheck::usage(usg);
+  // tmp output
+  p0 = inRange(0.00001, 0.9999, atof(usg.optionValue('p', 0).c_str()));
+  P  = inRange(0.00001, 0.9999, atof(usg.optionValue('P', 0).c_str()));
+  for (size_t i = 0, I = input_bitmap.size(); i < I; i++) {
+    for (size_t j = 0, J = input_bitmap.size(); j < J; j++) {
+      size_t num = 3u;
+      messagePrintf( 1, "ransac %lu --- %lu", i, j);
+      MyRansacCheck chk(&(fpsv[i]), &(fpsv[j]));
+      FeaturePointIndexPairs ret = ransac(pairs[i][j], chk, num, p0, P);
+      if (!ret.empty()) {
+        chk.rememberVCalc(ret);
+        FeaturePointIndexPairs ok(ret);
+        for (size_t k = 0, K = pairs[i][j].size(); k < K; k++) {
+          bool chk_again = true;
+          for (size_t l = 0, L = ret.size(); chk_again && l < L; l++) {
+            if (ret[l] == pairs[i][j][k]) {
+              chk_again = false;
+            }
+          }
+          if (chk_again && chk.ok(pairs[i][j][k])) {
+            ok.push_back(pairs[i][j][k]);
+          }
+        }
+        if (ok.size() >= num) pairs[i][j] = ok;
+        else                  pairs[i][j].clear();
+        messagePrintf(-1, "ok(%lu)", pairs[i][j].size());
+      }
+      else {
+        pairs[i][j].clear();
+        messagePrintf(-1, "empty");
+      }
+    }
+  }
+  messagePrintf(-1, "ok");
+  return true;
+}
+
+///////////////////////// **# prob module #** ////////////////////////
+bool prob_mod() {
+  q = inRange(0.00001, 0.99999, atof(usg.optionValue('q', 0).c_str()));
+  r = inRange(0.00001, 0.99999, atof(usg.optionValue('r', 0).c_str()));
+  Q = inRange(0.00001, 0.99999, atof(usg.optionValue('Q', 0).c_str()));
+  double m_ni = log(q * (1 - r)) - log(r * (1 - q));
+  double c    = log(Q) - log(1 - Q);
+  double m_nf = log(1 - r) - log(1 - q);
+  messagePrintf(1, "run prob_mod, ni * %.7f > %.7f + nf * %.7f ???",
+                m_ni, c, m_nf);
+  for (size_t i = 0, I = input_bitmap.size(); i < I; i++) {
+    for (size_t j = 0, J = input_bitmap.size(); j < J; j++) {
+      if (pairs[i][j].empty()) {
+        continue;
+      }
+      double ni = pairs[i][j].size(), nf = 0;
+      MyRansacCheck chk(&(fpsv[i]), &(fpsv[j]));
+      chk.rememberVCalc(pairs[i][j]);
+      for (size_t k = 0, K = fpsv[i].size(); k < K; k++) {
+        Vector2D<double> to(chk.to(Vector2D<double>(fpsv[i][k](0),
+                                                    fpsv[i][k](1))));
+        if (0 <= to.x() && to.x() <= (double)input_bitmap[j].width() &&
+            0 <= to.y() && to.y() <= (double)input_bitmap[j].height()) {
+          nf++;
+        }
+      }
+      if (ni * m_ni > c + m_nf * nf) {
+        messagePrintf(0, "accept %lu --- %lu", i, j);
+        messagePrintf(0,
+                      "%.0f * %.3f = %.3f ?? %.3f = %.3f + %.3f * %.0f",
+                      ni, m_ni, ni * m_ni, c + m_nf * nf, c, m_nf, nf);
+        continue;
+      }
+      else {
+        pairs[i][j].clear();
+      }
+    }
+  }
+  messagePrintf(-1, "ok");
+  return true;
+}
+
+///////////////////// **# group them together #** ////////////////////
+bool group() {
+  messagePrintf(1, "group");
+  // union
+  DisjointSet dsj(input_bitmap.size());
+  for (size_t i = 0, I = input_bitmap.size(); i < I; i++) {
+    for (size_t j = 0; j < I; j++) {
+      if(pairs[i][j].empty()) continue;
+      dsj.merge(i, j);
+    }
+  }
+  std::vector<size_t> root;
+  for (size_t i = 0, I = input_bitmap.size(); i < I; i++) {
+    if (dsj.root(i) == i) {
+      root.push_back(i);
+    }
+  }
+  // split into groups
+  outputs.resize(root.size());
+  for (size_t i = 0, I = root.size(); i < I; i++) {
+    messagePrintf(1, "Group %d", i);
+    std::vector<size_t> ids;
+    for (size_t j = 0, J = input_bitmap.size(); j < J; j++) {
+      if (dsj.root(j) != root[i]) continue;
+      outputs[i].cameras.push_back(Camera<RGBf_Space>());
+      outputs[i].cameras[outputs[i].cameras.size() - 1].photo(
+        Photo<RGBf_Space>(input_bitmap[j])
+      );
+      ids.push_back(j);
+      messagePrintf(0, "camera %lu from bitmap %lu",
+                    outputs[i].cameras.size() - 1, j);
+    }
+    for (size_t j = 0, J = ids.size(); j < J; ++j) {
+      for (size_t k = 0; k < J; ++k) {
+        if (j == k) continue;
+        size_t i1 = ids[j], i2 = ids[k];
+        if (pairs[i1][i2].empty()) continue;
+        outputs[i].edges.push_back(OutputSet::Edge(j, k));
+        size_t index = outputs[i].edges.size() - 1;
+        for (size_t n = 0, N = pairs[i1][i2].size(); n < N; ++n) {
+          outputs[i].edges[index].v1.push_back(
+            fpsv[i1][pairs[i1][i2][n].from.second]
+          );
+          outputs[i].edges[index].v2.push_back(
+            fpsv[i1][pairs[i1][i2][n].to.second]
+          );
+        }
+        messagePrintf(0, "Edge %lu---%lu, size = %lu",
+                      i1, i2, outputs[i].edges[index].v1.size());
+      }
+    }
+    std::sort(outputs[i].edges.begin(), outputs[i].edges.end());
+    messagePrintf(-1, "");
+  }
+  messagePrintf(-1, "ok");
+  return true;
+}
+
+////////////////////// **# boundle adjustment #** ////////////////////
+bool boundle() {
+  stop = inRange(0.01, 100000.0, atof(usg.optionValue('s', 0).c_str()));
+  messagePrintf(1, "boundle adjustment");
+  for (size_t i = 0, I = outputs.size(); i < I; i++) {
+    int id = 0;
+    std::set<size_t> in;
+    size_t i1 = outputs[i].edges[0].i1;
+    size_t i2 = outputs[i].edges[0].i2;
+    for (size_t j = 0, J = outputs[i].edges[0].v1.size(); j < J; ++j) {
+      outputs[i].cameras[i1].fixedPoints2DGet().identityPointAdd(
+        id, outputs[i].edges[0].v1[j]);
+      outputs[i].cameras[i2].fixedPoints2DGet().identityPointAdd(
+        id, outputs[i].edges[0].v2[j]);
+      id++;
+    }
+    in.insert(i1);
+    in.insert(i2);
+    double r_lst = Camera<RGBf_Space>::boundleAdjustment2D(
+      &(outputs[i].cameras),
+      stop
+    );
+    for (size_t j = 1, J = outputs[i].edges.size(); j < J; ++j) {
+      size_t best;
+      for (best = 0; best < J; ++best) {
+        if (in.find(outputs[i].edges[best].i1) == in.end() &&
+            in.find(outputs[i].edges[best].i2) == in.end()) continue;
+        break;
+      }
+      i1 = outputs[i].edges[best].i1;
+      i2 = outputs[i].edges[best].i2;
+      for (size_t j = 0, J = outputs[i].edges[best].v1.size(); j < J; ++j) {
+        outputs[i].cameras[i1].fixedPoints2DGet().identityPointAdd(
+          id, outputs[i].edges[best].v1[j]);
+        outputs[i].cameras[i2].fixedPoints2DGet().identityPointAdd(
+          id, outputs[i].edges[best].v2[j]);
+        id++;
+      }
+      in.insert(i1);
+      in.insert(i2);
+      std::vector<Camera<RGBf_Space> > tmp(outputs[i].cameras);
+      double r = Camera<RGBf_Space>::boundleAdjustment2D(&tmp, stop);
+      if (r > r_lst * 1.5) continue;
+      outputs[i].cameras = tmp;
+    }
+  }
+  messagePrintf(-1, "ok");
+  return true;
+}
+
+bool expand() {
+  o_radius = inRange(100.0, 1000000.0, atof(usg.optionValue('O', 0).c_str()));
+  output_bitmap.resize(outputs.size());
+  for (size_t i = 0, I = outputs.size(); i < I; ++i) {
+    WatchBall<RGBf_Space> wb;
+    wb.cameras(outputs[i].cameras);
+    output_bitmap[i] = wb.expand(o_radius);
+  }
+  return true;
+}
+
+////////////////////// **# Write to output file #** //////////////////
+bool output() {
+  messagePrintf(1, "Write images");
+  for (size_t i = 0; i < output_bitmap.size(); i++) {
+    size_t width  = output_bitmap[i].width ();
+    size_t height = output_bitmap[i].height();
+    cv::Mat img(height, width, CV_8UC3);
+    for (size_t x = 0; x < width; x++) {
+      for (size_t y = 0; y < height; y++) {
+        RGBi_Space tmp;
+        colorTransformate(output_bitmap[i].pixel(y, x), &tmp);
+        img.at<cv::Vec3b>(y, x)[0] = tmp.b();
+        img.at<cv::Vec3b>(y, x)[1] = tmp.g();
+        img.at<cv::Vec3b>(y, x)[2] = tmp.r();
+      }
+    }
+    std::string output_name(usg.optionValue('o', 0)
+                            + (output_bitmap.size() > 1
+                               ? stringPrintf("%lu", i)
+                               : "")
+                            + ".jpg");
+    messagePrintf(1, "Write to file '%s'", output_name.c_str());
+    if (imwrite(output_name, img) == false) {
+      messagePrintf(-1, "opencv fail, ignore");
+    }
+    else {
+      messagePrintf(-1, "%lux%lu, ok", width, height);
+    }
+  }
+  messagePrintf(-1, "ok");
+  return true;
+}
+
+//*
+bool tmp_output() {
+  output_bitmap = input_bitmap;
+  for (size_t i = 0, I = input_bitmap.size(); i < I; i++) {
+    for (size_t j = 0, J = fpsv[i].size(); j < J; j++) {
+      ssize_t x = fpsv[i][j](0);
+      ssize_t y = fpsv[i][j](1);
+      ssize_t dx[2] = {0, 1}, x0[2] = {0, -10};
+      ssize_t dy[2] = {1, 0}, y0[2] = {-10, 0};
+      for(size_t k = 0; k < 2; k++){
+        for(size_t count = 0; count < 20; count++){
+          ssize_t xx = x + dx[k] * count + x0[k];
+          ssize_t yy = y + dy[k] * count + y0[k];
+          if(0 <= xx && xx < (ssize_t)input_bitmap[i].width() &&
+             0 <= yy && yy < (ssize_t)input_bitmap[i].height()){
+            output_bitmap[i].pixel(yy, xx, Vector3D<double>(1.0, 1.0, 0.0));
+          }
+        }
+      }
+    }
+  }
+  return output();
+}
+// */
+
+/*
+bool g_output(){
+  output_bitmap.resize(input_bitmap.size() * 2);
+  for(size_t i = 0, I = input_bitmap.size(); i < I; i++){
+    output_bitmap[i * 2    ] = input_bitmap[i];
+    output_bitmap[i * 2 + 1] = input_bitmap[i];
+    output_bitmap[i * 2    ].gradiancedX(3, 3);
+    output_bitmap[i * 2 + 1].gradiancedY(3, 3);
+    for(size_t x = 0, X = output_bitmap[i * 2].width(); x < X; x++){
+      for(size_t y = 0, Y = output_bitmap[i * 2].height(); y < Y; y++){
+        Vector3D<double> v;
+        v = output_bitmap[i * 2](y, x);
+        output_bitmap[i * 2](y, x) = Vector3D<double>(v.length() / sqrt(2.0));
+        v = output_bitmap[i * 2 + 1](y, x);
+        output_bitmap[i * 2 + 1](y, x) = Vector3D<double>(v.length() / sqrt(2.0));
+      }
+    }
+  }
+  return output();
+}
+
+// */
+
+bool pair_output(){
+  for(size_t i = 0, I = input_bitmap.size(); i < I; i++){
+    for(size_t j = 0, J = input_bitmap.size(); j < J; j++){
+      if(pairs[i][j].empty()) continue;
+      MyRansacCheck chk(&(fpsv[i]), &(fpsv[j]));
+      chk.rememberVCalc(pairs[i][j]);
+      size_t index = output_bitmap.size();
+      output_bitmap.push_back(input_bitmap[i]);
+      for(ssize_t x = 0, X = input_bitmap[i].width(); x < X; x++){
+        for(ssize_t y = 0, Y = input_bitmap[i].height(); y < Y; y++){
+          Vector2D<double> to(chk.to(Vector2D<double>(x, y)));
+          ssize_t x2 = to.x(), y2 = to.y();
+          if(0 <= x2 && x2 <= (ssize_t)input_bitmap[j].width() &&
+             0 <= y2 && y2 <= (ssize_t)input_bitmap[j].height()){
+            output_bitmap[index].pixel(y, x, (input_bitmap[i].pixel(y, x) +
+                                              input_bitmap[j].pixel(y2,x2)) / 2
+                                       );
+          }
+        }
+      }
+    }
+  }
+  return output();
+}
+/*
+bool pair_output2(){
+  for(size_t i = 0, I = input_bitmap.size(); i < I; i++){
+    for(size_t j = 0, J = input_bitmap.size(); j < J; j++){
+      if((i + 1) % I != j && (j + 1) % J != i) continue;
+      messagePrintf(0, "%3lu--%3lu: %lu", i, j, pairs[i][j].size());
+      if(pairs[i][j].empty()) continue;
+      MyRansacCheck chk(&(fpsv[i]), &(fpsv[j]));
+      chk.rememberVCalc(pairs[i][j]);
+      size_t index = output_bitmap.size();
+      output_bitmap.push_back(input_bitmap[i]);
+      for(ssize_t x = 0, X = input_bitmap[i].width(); x < X; x++){
+        for(ssize_t y = 0, Y = input_bitmap[i].height(); y < Y; y++){
+          Vector2D<double> to(chk.to(Vector2D<double>(x, y)));
+          ssize_t x2 = to.x(), y2 = to.y();
+          if(0 <= x2 && x2 <= (ssize_t)input_bitmap[j].width() &&
+             0 <= y2 && y2 <= (ssize_t)input_bitmap[j].height()){
+            output_bitmap[index].pixel(y, x, (input_bitmap[i].pixel(y, x) +
+                                              input_bitmap[j].pixel(y2,x2)) / 2
+                                       );
+          }
+        }
+      }
+      for(size_t k = 0, K = fpsv[i].size(); k < K; k++){
+        ssize_t dy[2] = {0, 1}, dx[2] = {1, 0};
+        ssize_t y0[2] = {0, -10}, x0[2] = {-10, 0};
+        ssize_t x = fpsv[i][k](0), y = fpsv[i][k](1);
+        for(ssize_t m = 0; m < 2; m++){
+          for(ssize_t n = 0; n < 20; n++){
+            ssize_t xx = x + x0[m] + dx[m] * n;
+            ssize_t yy = y + y0[m] + dy[m] * n;
+            if(0 <= xx && xx <= (ssize_t)input_bitmap[j].width() &&
+               0 <= yy && yy <= (ssize_t)input_bitmap[j].height()){
+              output_bitmap[index].pixel(yy, xx,
+                                         Vector3D<double>(1.0, 1.0, 0));
+            }
+          }
+        }
+      }
+      for(size_t k = 0, K = pairs[i][j].size(); k < K; k++){
+        ssize_t dy[2] = {0, 1}, dx[2] = {1, 0};
+        ssize_t y0[2] = {0, -10}, x0[2] = {-10, 0};
+        ssize_t x = fpsv[i][pairs[i][j][k].from.second](0);
+        ssize_t y = fpsv[i][pairs[i][j][k].from.second](1);
+        for(ssize_t m = 0; m < 2; m++){
+          for(ssize_t n = 0; n < 20; n++){
+            ssize_t xx = x + x0[m] + dx[m] * n;
+            ssize_t yy = y + y0[m] + dy[m] * n;
+            if(0 <= xx && xx <= (ssize_t)input_bitmap[j].width() &&
+               0 <= yy && yy <= (ssize_t)input_bitmap[j].height()){
+              output_bitmap[index].pixel(yy, xx,
+                                         Vector3D<double>(1.0, 0.0, 0));
+            }
+          }
+        }
+      }
+    }
+  }
+  return output();
+}
+// */
+
+
+int main(int argc, char** argv){
+  setup(argc, argv);
+  input();
+  detect();
+  kmatch();
+  ransac();
+  prob_mod();
+  group();
+  pair_output(); return 0;
+  boundle();
+  expand();
+  output();
+  return 0;
+}