ExampleRANSAC_double.cpp

/* 
This file is part of the BIAS library (Basic ImageAlgorithmS).

Copyright (C) 2003-2009    (see file CONTACT for details)
  Multimediale Systeme der Informationsverarbeitung
  Institut fuer Informatik
  Christian-Albrechts-Universitaet Kiel


BIAS is free software; you can redistribute it and/or modify
it under the terms of the GNU Lesser General Public License as published by
the Free Software Foundation; either version 2.1 of the License, or
(at your option) any later version.

BIAS is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
GNU Lesser General Public License for more details.

You should have received a copy of the GNU Lesser General Public License
along with BIAS; if not, write to the Free Software
Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
*/

/** @example ExampleRANSAC_double.cpp
    @relates DoubleRANSAC, RANSAC
    @ingroup g_examples    
    @brief   Example for using classes derived from BIAS::RANSAC to
             solve custom problems robustly.
    @author MIP
*/

#include <Base/Common/BIASpragma.hh>

#include "ExampleRANSAC_double.hh"

#include <Base/Math/Random.hh>
#include <vector>
#include <iostream>
#include <iomanip>
#include <float.h>

using namespace BIAS;
using namespace std;

/* Generate test data for RANSAC_double instance.
   input : inlier   percentage of inliers
           count    number of data generated
           mean     mean ofthe generated data
           sigma    standard deviation of the generated data
   output: data     the generated data
*/
void GenerateRANSAC_doubleData(unsigned int count, double mean, double sigma, 
                               double inlier_fraction, vector<double>& data, 
                               Random& randomizer)
{
  if (data.size() != count) data.resize(count);
  double outlier_prob;
  for (unsigned int i = 0; i < count; i++) {
    outlier_prob = randomizer.GetUniformDistributed(0.0, 1.0);
    if (outlier_prob > inlier_fraction) {
      // generate a outlier
      data[i] = randomizer.GetNormalDistributed(mean, sigma) +
                randomizer.GetUniformDistributed(3*sigma, 30.0*sigma);
    } else {
      // generate an inlier
      data[i] = randomizer.GetNormalDistributed(mean, sigma);
    }
  }
}


int main(int argc, char *argv[])
{
  Random randomizer;
  int res = 0;
  vector<vector<double> > data_array;
  vector<double> data;
  double min = -100.0, max = 100.0;
  double minsigma = 0.0, maxsigma = 5.0;
  vector<double> mean;
  vector<double> sigma;
  double inlier_fraction = 0.7;
  vector<double> expected_inlier_fraction;
  unsigned int sample_size = 2;
  bool refine_solution = true;
  double inlier_distance_threshold;
  double solution;
  vector<bool> inliers;
  bool debug = false; //true;
  RANSAC_double ransac;
  double mean_error = 0.0;
  double solution_error_threshold;
  
  if (debug) {
    ransac.AddDebugLevel("D_RANSAC_SAMPLES");
    ransac.AddDebugLevel("D_RANSAC_SOLVE_MASTER");
    ransac.AddDebugLevel("D_DOUBLE_RANSAC_SOLUTION");
  }

  unsigned int count = argc > 1 ? atoi(argv[1]) : 1000;
  unsigned int data_length = argc > 2 ? atoi(argv[2]) : 1000;
  cout << endl << "Running " << argv[0] << " " << count << " times with "
       << data_length << " samples per test : " << endl
       << "(call " << argv[0] << " <num_tests> <num_samples> to change "
       << "parameters)" << endl << endl;
  
  // generate all test sets in advance since randomizer might be reset
  // in RANSAC instance
  data_array.resize(count);
  mean.resize(count);
  sigma.resize(count);
  expected_inlier_fraction.resize(count);
  for (unsigned int i = 0; i < count; i++)
  {
    mean[i] = randomizer.GetUniformDistributed(min, max);
    sigma[i] = randomizer.GetUniformDistributed(minsigma, maxsigma);
    expected_inlier_fraction[i] = inlier_fraction +
      randomizer.GetNormalDistributed(0.0, (1.0-inlier_fraction)/3.0);
    if (expected_inlier_fraction[i] > 1.0)
      expected_inlier_fraction[i] = 1.0;
    GenerateRANSAC_doubleData(data_length, mean[i], sigma[i], inlier_fraction, 
                              data_array[i], randomizer);
  }

  // evaluate all test sets
  for (unsigned int i = 0; i < count; i++)
  {
    inlier_distance_threshold = 3.0*sigma[i];
    solution_error_threshold = 1.5*sigma[i];
    //expected_inlier_fraction[i] = 0.5;
    if (debug) {
      for (unsigned int j = 0; j < data_length; j++)
        cout << data_array[i][j] << "  ";
      cout << endl;
      cout <<setw(5)<<i<<":\tmean: "<<mean[i]<<"\tsigma: "<<sigma[i]<<endl;
      cout << "inlier_distance_threshold = "<<inlier_distance_threshold<<endl;
      cout << "solution_error_threshold = "<<solution_error_threshold<<endl;
      cout << "expected_inlier_fraction " << expected_inlier_fraction[i]<<endl;
    }
  
    ransac.Init(data_array[i], sample_size, inlier_distance_threshold, 
                refine_solution, expected_inlier_fraction[i], 
                solution_error_threshold);

    if (ransac.SolveMaster(expected_inlier_fraction[i], solution, inliers) < 0)
    {
      BIASERR("Error running RANSAC_double::SolveMaster()!");
      cerr << "- true solution          \t" << mean[i] <<"  \t"
           << "sigma: "<<sigma[i]<<"    \t"
           << "expected_inlier_fraction: "<<expected_inlier_fraction[i]<<endl;
    } else {
      cerr << "- true solution          \t" << mean[i] <<"  \t"
           << "sigma: "<<sigma[i]<<"    \t"
           << "expected_inlier_fraction: "<<expected_inlier_fraction[i]<<endl
           << "  estimated solution     \t" << solution <<"  \t"
           << "error: "<<fabs(solution-mean[i])<<"    \t"
           << "expected_inlier_fraction: "<<expected_inlier_fraction[i]<<endl;
      mean_error += (sigma[i]!=0.0) ? (fabs(solution-mean[i])/sigma[i]) :
                    (fabs(solution-mean[i]) / DBL_MIN);
    }
  }
  mean_error /= (double)count;

  cerr << "\n mean error of all tests : " << mean_error << endl << endl;
  return res;
}