Commit a2ee7d45 authored by Laurent Belcour's avatar Laurent Belcour

Corrections on the rational 1D fitter + adding parsing of the arguments

parent 3afec497
#pragma once
#include <string>
#include <map>
#include <cstdlib>
#include <iostream>
class arguments
{
public: // functions
// Constructor and destructor
arguments(int argc, char** const argv)
{
std::string key ;
std::string data ;
for(int i=0; i<argc; ++i)
{
std::string temp(argv[i]) ;
std::string data ;
if(temp.compare(0, 2, "--") == 0)
{
key = temp.substr(2, temp.size()-2) ;
if(i+1 < argc)
{
std::string next(argv[i+1]) ;
if(next.compare(0, 2, "--") != 0)
{
data = next ;
}
}
}
_map.insert(std::pair<std::string, std::string>(key, data)) ;
}
} ;
~arguments() { } ;
// Access elements
bool is_defined(const std::string& key) const
{
if(_map.count(key) > 0)
{
return true ;
}
else
{
return false ;
}
} ;
std::string operator[](const std::string& key) const
{
if(_map.count(key) > 0)
{
return _map.find(key)->second ;
}
else
{
std::cerr << "Underfined request to key : \"" << key << "\"" << std::endl ;
return std::string() ;
}
} ;
float get_float(const std::string& key, float default_value = 0.0f)
{
std::string value = _map[key] ;
if(value.empty())
return default_value ;
else
return atof(value.c_str()) ;
} ;
int get_int(const std::string& key, int default_value = 0)
{
std::string value = _map[key] ;
if(value.empty())
return default_value ;
else
return atoi(value.c_str()) ;
} ;
private: // data
std::map<std::string, std::string> _map ;
} ;
......@@ -4,8 +4,13 @@
#include <CGAL/QP_models.h>
#include <CGAL/QP_functions.h>
#include <CGAL/MP_Float.h>
#include <boost/regex.hpp>
#include <string>
#include <iostream>
#include <fstream>
#include <limits>
#include <algorithm>
rational_1d::rational_1d()
{
......@@ -57,44 +62,117 @@ void rational_1d::save() const
{
}
std::ostream& operator<< (std::ostream& out, const rational_1d& r)
{
std::cout << "p = [" ;
for(int i=0; i<r.a.size(); ++i)
{
if(i != 0)
{
std::cout << ", " ;
}
std::cout << r.a[i] ;
}
std::cout << "]" << std::endl ;
std::cout << "q = [" ;
for(int i=0; i<r.b.size(); ++i)
{
if(i != 0)
{
std::cout << ", " ;
}
std::cout << r.b[i] ;
}
std::cout << "]" << std::endl ;
}
void rational_1d_data::load(const std::string& filename)
{
std::ifstream file(filename) ;
_min = std::numeric_limits<float>::max() ;
_max = -std::numeric_limits<float>::max() ;
if(!file.is_open())
{
std::cerr << "<<ERROR>> unable to open file \"" << filename << "\"" << std::endl ;
}
// N-Floats regexp
boost::regex e ("^([0-9]*\.?[0-9]+[\\t ]?)+");
float x, y, dy ;
while(file.good())
{
file >> x >> y >> dy ;
std::string line ;
std::getline(file, line) ;
std::stringstream linestream(line) ;
// Discard incorrect lines
if(!boost::regex_match(line,e))
{
continue ;
}
linestream >> x >> y ;
if(linestream.good()) {
linestream >> dy ;
} else {
// TODO Specify the delta in case
dy = 0.1f ;
}
std::vector<float> v ;
v.push_back(x) ;
v.push_back(y-dy) ;
v.push_back(y+dy) ;
data.push_back(v) ;
_data.push_back(v) ;
// Update min and max
_min = std::min(_min, x) ;
_max = std::max(_max, x) ;
}
// Sort the vector
std::sort(_data.begin(), _data.end(), [](const std::vector<float>& a, const std::vector<float>& b){return (a[0]<b[0]);});
std::cout << "<<INFO>> loaded file \"" << filename << "\"" << std::endl ;
std::cout << "<<INFO>> data inside [" << _min << ", " << _max << "]" << std::endl ;
}
bool rational_1d_data::get(int i, float& x, float& yl, float& yu) const
{
if(i >= data.size())
if(i >= (int)_data.size())
{
return false ;
}
x = data[i][0] ;
yl = data[i][1] ;
yu = data[i][2] ;
x = _data[i][0] ;
yl = _data[i][1] ;
yu = _data[i][2] ;
return true ;
}
const std::vector<float>& rational_1d_data::operator[](int i) const
{
return data[i] ;
return _data[i] ;
}
int rational_1d_data::size() const
{
return data.size() ;
return _data.size() ;
}
float rational_1d_data::min() const
{
return _min ;
}
float rational_1d_data::max() const
{
return _max ;
}
typedef CGAL::MP_Float ET ;
......@@ -103,6 +181,11 @@ typedef CGAL::Quadratic_program_solution<ET> Solution ;
// Fitting a data
rational_1d rational_1d_fitter::fit_data(const rational_1d_data& data)
{
return fit_data(data, 10, 10) ;
}
rational_1d rational_1d_fitter::fit_data(const rational_1d_data& data, int np, int nq)
{
// The resulting ration function
rational_1d r ;
......@@ -112,12 +195,11 @@ rational_1d rational_1d_fitter::fit_data(const rational_1d_data& data)
// Select the size of the result vector to
// be equal to the dimension of p + q
int np = 10, nq = 10;
for(int i=0; i<np+nq; ++i)
{
qp.set_d(i, i, 1.0f) ;
}
// Each constraint (fitting interval or point
// add another dimension to the constraint
// matrix
......@@ -168,12 +250,13 @@ rational_1d rational_1d_fitter::fit_data(const rational_1d_data& data)
if(qp.is_linear() && !qp.is_nonnegative())
{
std::cerr << "Error here, should not be linear or negative!" << std::endl ;
std::cerr << "<<ERROR>> the problem should not be linear or negative!" << std::endl ;
throw ;
}
#ifdef DEBUG
#ifdef EXTREM_DEBUG
// Iterate over the rows
std::cout << std::endl ;
std::cout << "A = [" ;
for(int j=0; j<qp.get_m(); ++j)
{
......@@ -188,13 +271,33 @@ rational_1d rational_1d_fitter::fit_data(const rational_1d_data& data)
std::cout << ";" ;
if(j < qp.get_n()-1) std::cout << std::endl ;
}
std::cout << "]" << std::endl ;
std::cout << "]" << std::endl << std::endl ;
std::cout << std::endl ;
std::cout << "D = [" ;
for(int j=0; j<np+nq; ++j)
{
if(j == 0) std::cout << " " ;
// Iterate over the columns
for(int i=0; i<np+nq; ++i)
{
if(i > 0) std::cout << ",\t" ;
std::cout << *(*(qp.get_d()+i) +j) ;
}
std::cout << ";" ;
}
std::cout << "]" << std::endl << std::endl ;
#endif
// solve the program, using ET as the exact type
//*
Solution s = CGAL::solve_quadratic_program(qp, ET()) ;
assert (s.solves_quadratic_program(qp)) ;
/*/
Solution s = CGAL::solve_nonnegative_quadratic_program(qp, ET());
assert (s.solves_nonnegative_quadratic_program(qp));
//*/
// Recopy the vector data
std::vector<float> p, q;
for(int i=0; i<np+nq; ++i)
......
......@@ -28,6 +28,9 @@ class rational_1d : public std::function<float(float)>
void load(const std::string& filename) ;
void save() const ;
// STL stream ouput
friend std::ostream& operator<< (std::ostream& out, const rational_1d& r) ;
private: // data
// Store the coefficients for the moment, I assume
......@@ -50,11 +53,19 @@ class rational_1d_data // : public fitting_data
// Get data size
int size() const ;
// Get min and max input parameters
float min() const ;
float max() const ;
private: // data
// Store for each point of data, the upper
// and lower value
std::vector<std::vector<float> > data ;
std::vector<std::vector<float> > _data ;
// Store the min and max value on the input
// domain
float _min, _max ;
} ;
class rational_1d_fitter // : public fitting_algorithm
......@@ -64,5 +75,9 @@ class rational_1d_fitter // : public fitting_algorithm
// Fitting a data object
rational_1d fit_data(const rational_1d_data& data) ;
// Fitting a data object using np elements
// in the numerator and nq elements in the
// denominator
rational_1d fit_data(const rational_1d_data& data, int np, int nq) ;
} ;
......@@ -8,8 +8,8 @@ int main(int argc, int argv)
for(int i=0; i<100; ++i)
{
const float x = i / (float)99.0f ;
const float y = 100.0f * exp(-100.0f * x*x) * x*x ;
const float x = i / (float)10.0f ;
const float y = exp(-0.2 * x*x) * x*x + 0.01 * (sin(0.3*x) * x*x);
f << x << "\t" << y << "\t" << 0.01f << std::endl ;
}
......
......@@ -4,17 +4,36 @@
#include <iostream>
#include <fstream>
#include <core/args.h>
int main(int argc, char** argv)
{
arguments args(argc, argv) ;
if(args.is_defined("help")) {
std::cout << argv[0] << " --np <int> --nq <int> --input <filename> --output <filename>" << std::endl ;
return 0 ;
}
if(! args.is_defined("input")) {
std::cerr << "<<ERROR>> the input filename is not defined" << std::endl ;
return 1 ;
}
if(! args.is_defined("output")) {
std::cerr << "<<ERROR>> the output filename is not defined" << std::endl ;
return 1 ;
}
rational_1d_data data ;
data.load(std::string("input.gnuplot"));
data.load(args["input"]);
rational_1d_fitter fitter ;
rational_1d r = fitter.fit_data(data) ;
rational_1d r = fitter.fit_data(data, args.get_int("np", 10), args.get_int("nq", 10)) ;
std::cout << r << std::endl ;
//*
std::ofstream file("output.gnuplot", std::ios_base::trunc);
for(float x=0.0f; x<=1.09f; x+=0.01f)
std::ofstream file(args["output"], std::ios_base::trunc);
const float dt = (data.max() - data.min()) / 100.0f ;
for(float x=data.min(); x<=data.max(); x+=dt)
{
file << x << "\t" << r(x) << std::endl ;
}
......
DEST_DIR = ../bin
CONFIG += debug
INCLUDEPATH += ../../libs/rational_1d
INCLUDEPATH += ../../ ../../libs/rational_1d
SOURCES += main.cpp ../../libs/rational_1d/rational_1d_fitter.cpp
QMAKE_CXXFLAGS += -std=c++11 -frounding-math
LIBS += -lCGAL
LIBS += -lCGAL -lboost_regex
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment