D:/Programmierung/C/TRTK/include/TRTK/Ransac.hpp

Go to the documentation of this file.

/*
    Copyright (C) 2010 - 2014 Christoph Haenisch

    Chair of Medical Engineering (mediTEC)
    RWTH Aachen University
    Pauwelsstr. 20
    52074 Aachen
    Germany

    See license.txt for more information.

    Version 0.1.0 (2012-03-23)
*/

#ifndef RANSAC_HPP_1147104378
#define RANSAC_HPP_1147104378


#include <cstddef>
#include <cmath>
#include <limits>
#include <vector>
#include <set>

#include "Coordinate.hpp"
#include "ErrorObj.hpp"


namespace TRTK
{


template <class T, class DataType = Coordinate<T> >
class Ransac
{
public:

    class Model
    {
    public:
        virtual void compute() = 0;                                     
        virtual T getDeviation(const DataType & datum) const = 0;       
        virtual unsigned getMinimumNumberOfItems() const = 0;           
        virtual T getRMSE() const = 0;                                  
        virtual void setData(const std::vector<DataType> & data) = 0;   
    };                                                                  

    enum Algorithm
    {
        RANSAC,                             
        RANSAC_BIGGEST_CONSESUS_SET,        
        RANSAC_SMALLEST_RMSE                
    };

    enum Error
    {
        NO_DATA_AVAILABLE,                  
        NO_MODEL_AVAILABLE,                 
        NOT_ENOUGH_DATA_POINTS,             
        UNKNOWN_ALGORITHM,                  
        UNKNOWN_ERROR                       
    };

    Ransac();
    virtual ~Ransac();

    unsigned compute();

    void setAlgorithm(Algorithm algorithm = RANSAC);
    void setData(const std::vector<DataType> & data);
    void setErrorTolerance(T error_tolerance);
    void setMinimumSetSize(unsigned threshold);
    void setModel(Model & model);
    void setProbability(T probability);
    void setThreshold(unsigned threshold);
    void setTrials(unsigned number);

protected:
    const std::vector<DataType> * data;
    Model * model;

    T error_tolerance;
    unsigned threshold;
    unsigned trials;
    T probability;                          // Probability that any selected data point of a minimal subset of data is within the error tolerance of the model.
    unsigned minimumSetSize;

    unsigned (Ransac::*algorithm)();

    unsigned algorithmRansac();
    unsigned algorithmRansacBiggestConsensusSet();
    unsigned algorithmRansacSmallestRmse();
};


template <class T, class DataType>
Ransac<T, DataType>::Ransac() :
    data(NULL),
    model(NULL),
    error_tolerance(0),
    threshold(0),
    trials(0),
    probability(0.2),
    minimumSetSize(0),
    algorithm(&Ransac::algorithmRansac)
{
}


template <class T, class DataType>
Ransac<T, DataType>::~Ransac()
{
}


template <class T, class DataType>
unsigned Ransac<T, DataType>::algorithmRansac()
{
    using namespace std;
    using namespace Tools;

    const unsigned maximum_number_of_items = data->size();
    const unsigned minimum_number_of_items = model->getMinimumNumberOfItems();

    vector<DataType> initial_data_set(minimum_number_of_items);
    vector<DataType> consensus_set;
    vector<DataType> best_consensus_set;

    // For each trial...

    for (unsigned i = 0; i < trials; ++i)
    {
        // Randomly select items from the data set and use them
        // as an initial data set for the model fitting.

        set<int> indices;

        for (unsigned j = 0; j < minimum_number_of_items; ++j)
        {
            // Avoid using the same data point several times.

            unsigned index = rand(0u, maximum_number_of_items - 1);;

            while (indices.find(index) != indices.end())
            {
                index = rand(0u, maximum_number_of_items - 1);
            }

            indices.insert(index);

            // Store the randomly selected data point.

            initial_data_set[j] = (*data)[index];
        }

        // Estimate a model from the initial data set.

        model->setData(initial_data_set);
        model->compute();

        // Construct a new consensus set and add those items from the
        // original data set which lie within some error tolerance.

        consensus_set.clear();
        consensus_set.reserve(maximum_number_of_items);

        for (unsigned j = 0; j < maximum_number_of_items; ++j)
        {
            if (model->getDeviation((*data)[j]) < error_tolerance)
            {
                consensus_set.push_back((*data)[j]);
            }
        }

        if (consensus_set.size() <= minimum_number_of_items)
        {
            // Actually, this should never happen!
            // But you never know...

            consensus_set = initial_data_set;
        }

        // If the size of the consensus set is greater than the given
        // threshold, then compute a new model from this set and exit.

        if (consensus_set.size() > threshold)
        {
            model->setData(consensus_set);
            model->compute();

            return consensus_set.size();
        }

        // The size of the consensus set is *not* greater than the given
        // threshold, so randomly select a new subset. Save the old
        // consensus set in case no better one is found.

        if (best_consensus_set.size() < consensus_set.size())
        {
            best_consensus_set = consensus_set;
        }
    }

    // Use the consensus set to compute a new model.

    model->setData(best_consensus_set);
    model->compute();

    return best_consensus_set.size();
}


template <class T, class DataType>
unsigned Ransac<T, DataType>::algorithmRansacBiggestConsensusSet()
{
    using namespace std;
    using namespace Tools;

    const unsigned maximum_number_of_items = data->size();
    const unsigned minimum_number_of_items = model->getMinimumNumberOfItems();

    vector<DataType> initial_data_set(minimum_number_of_items);
    vector<DataType> consensus_set;
    vector<DataType> best_consensus_set;

    // For each trial...

    for (unsigned i = 0; i < trials; ++i)
    {
        // Randomly select items from the data set and use them
        // as an initial data set for the model fitting.

        set<int> indices;

        for (unsigned j = 0; j < minimum_number_of_items; ++j)
        {
            // Avoid using the same data point several times.

            unsigned index = rand(0u, maximum_number_of_items - 1);;

            while (indices.find(index) != indices.end())
            {
                index = rand(0u, maximum_number_of_items - 1);
            }

            indices.insert(index);

            // Store the randomly selected data point.

            initial_data_set[j] = (*data)[index];
        }

        // Estimate a model from the initial data set.

        model->setData(initial_data_set);
        model->compute();

        // Construct a new consensus set and add those items from the
        // original data set which lie within some error tolerance.

        consensus_set.clear();
        consensus_set.reserve(maximum_number_of_items);

        for (unsigned j = 0; j < maximum_number_of_items; ++j)
        {
            if (model->getDeviation((*data)[j]) < error_tolerance)
            {
                consensus_set.push_back((*data)[j]);
            }
        }

        if (consensus_set.size() <= minimum_number_of_items)
        {
            // Actually, this should never happen!
            // But you never know...

            consensus_set = initial_data_set;
        }

        // Randomly select a new subset. Save the old consensus set
        // in case no better one is found.

        if (best_consensus_set.size() < consensus_set.size())
        {
            best_consensus_set = consensus_set;
        }
    }

    // Use the biggest consensus set to compute a new model.

    model->setData(best_consensus_set);
    model->compute();

    return best_consensus_set.size();
}


template <class T, class DataType>
unsigned Ransac<T, DataType>::algorithmRansacSmallestRmse()
{
    using namespace std;
    using namespace Tools;

    const unsigned maximum_number_of_items = data->size();
    const unsigned minimum_number_of_items = model->getMinimumNumberOfItems();

    vector<DataType> initial_data_set(minimum_number_of_items);
    vector<DataType> consensus_set;
    vector<DataType> best_consensus_set;

    T rmse = numeric_limits<T>::max();

    // For each trial...

    for (unsigned i = 0; i < trials; ++i)
    {
        // Randomly select items from the data set and use them
        // as an initial data set for the model fitting.

        set<int> indices;

        for (unsigned j = 0; j < minimum_number_of_items; ++j)
        {
            // Avoid using the same data point several times.

            unsigned index = rand(0u, maximum_number_of_items - 1);;

            while (indices.find(index) != indices.end())
            {
                index = rand(0u, maximum_number_of_items - 1);
            }

            indices.insert(index);

            // Store the randomly selected data point.

            initial_data_set[j] = (*data)[index];
        }

        // Estimate a model from the initial data set.

        model->setData(initial_data_set);
        model->compute();

        // Construct a new consensus set and add those items from the
        // original data set which lie within some error tolerance.

        consensus_set.clear();
        consensus_set.reserve(maximum_number_of_items);

        for (unsigned j = 0; j < maximum_number_of_items; ++j)
        {
            if (model->getDeviation((*data)[j]) < error_tolerance)
            {
                consensus_set.push_back((*data)[j]);
            }
        }

        if (consensus_set.size() <= minimum_number_of_items)
        {
            // Actually, this should never happen!
            // But you never know...

            consensus_set = initial_data_set;
        }

        // Compute the model with the current consensus set.

        model->setData(consensus_set);
        model->compute();

        T rmse_current_model = model->getRMSE();

        // Randomly select a new subset. Save the old consensus set
        // in case no better one is found.

        if (rmse_current_model < rmse && consensus_set.size() >= minimumSetSize)
        {
            best_consensus_set = consensus_set;
        }
    }

    // Use the biggest consensus set to compute a new model.

    if (best_consensus_set.size() > minimum_number_of_items)
    {
        model->setData(best_consensus_set);
        model->compute();

        return best_consensus_set.size();
    }
    else
    {
        return 0;
    }
}


template <class T, class DataType>
unsigned Ransac<T, DataType>::compute()
{
    // Checking some assumptions.

    if (!model)
    {
        ErrorObj error;
        error.setClassName("Ransac");
        error.setFunctionName("compute");
        error.setErrorMessage("Regression model was not set.");
        error.setErrorCode(NO_MODEL_AVAILABLE);
        throw error;
    }

    if (!data)
    {
        ErrorObj error;
        error.setClassName("Ransac");
        error.setFunctionName("compute");
        error.setErrorMessage("No data available.");
        error.setErrorCode(NO_DATA_AVAILABLE);
        throw error;
    }

    if (data->size() < model->getMinimumNumberOfItems())
    {
        ErrorObj error;
        error.setClassName("Ransac");
        error.setFunctionName("compute");
        error.setErrorMessage("Not enough data points to compute the model.");
        error.setErrorCode(NOT_ENOUGH_DATA_POINTS);
        throw error;
    }

    using namespace std;

    // Compute the maximum number of attempts (trials) to find a consensus
    // set if not given. See the paper of Fischler and Bolles.

    unsigned trials_backup = trials;

    if (trials == 0)
    {
        int n = model->getMinimumNumberOfItems();
        trials = unsigned(ceil(3.0 * pow(probability, -T(n))));
    }

    // If no threshold is given, set it to 8/10 of the size of "data".

    unsigned threshold_backup = threshold;

    if (threshold == 0)
    {
        threshold = unsigned(ceil(0.8 * data->size()));
    }

    // If the minimum set size is greater than the size of "data" shrink
    // it to #data - model->getMinimumNumberOfItems().

    unsigned minimumSetSize_backup = minimumSetSize;

    if (minimumSetSize > data->size() && data->size() > model->getMinimumNumberOfItems())
    {
        minimumSetSize = data->size() - model->getMinimumNumberOfItems();
    }

    // Start the computation.

    unsigned number_of_items_used = (this->*algorithm)();

    // Restore the changed variables.

    minimumSetSize = minimumSetSize_backup;
    threshold = threshold_backup;
    trials = trials_backup;

    return number_of_items_used;
}


template <class T, class DataType>
void Ransac<T, DataType>::setAlgorithm(Algorithm algorithm)
{
    switch(algorithm)
    {
        case RANSAC:
        {
            this->algorithm = &Ransac::algorithmRansac;
            break;
        }

        case RANSAC_BIGGEST_CONSESUS_SET:
        {
            this->algorithm = &Ransac::algorithmRansacBiggestConsensusSet;
            break;
        }

        case RANSAC_SMALLEST_RMSE:
        {
            this->algorithm = &Ransac::algorithmRansacSmallestRmse;
            break;
        }

        default:
        {
            ErrorObj error;
            error.setClassName("Ransac");
            error.setFunctionName("setAlgorithm");
            error.setErrorMessage("Unknown algorithm.");
            error.setErrorCode(UNKNOWN_ALGORITHM);
            throw error;
        }
    }
}


template <class T, class DataType>
void Ransac<T, DataType>::setData(const std::vector<DataType> & data)
{
    this->data = &data;
}


template <class T, class DataType>
void Ransac<T, DataType>::setErrorTolerance(T error_tolerance)
{
    this->error_tolerance = error_tolerance;
}


template <class T, class DataType>
void Ransac<T, DataType>::setMinimumSetSize(unsigned size)
{
    minimumSetSize = size;
}


template <class T, class DataType>
void Ransac<T, DataType>::setModel(Model & model)
{
    this->model = &model;
}


template <class T, class DataType>
void Ransac<T, DataType>::setProbability(T probability)
{
    this->probability = probability;
}


template <class T, class DataType>
void Ransac<T, DataType>::setThreshold(unsigned threshold)
{
    this->threshold = threshold;
}


template <class T, class DataType>
void Ransac<T, DataType>::setTrials(unsigned number)
{
    trials = number;
}


} // namespace TRTK


#endif // RANSAC_HPP_1147104378