Commit 9ae66651 authored by Laurent Belcour's avatar Laurent Belcour

First attempt to add a non linear fitter using the CERES library from Google inc.

parent db788d74
......@@ -427,7 +427,7 @@ class compound_function: public nonlinear_function, public std::vector<nonlinear
virtual void save_call(std::ostream& out, const arguments& args) const
{
bool is_cpp = args["export"] == "C++";
bool is_shader = args["export"] == "shader";
bool is_shader = args["export"] == "shader" || args["export"] == "explorer";
bool is_matlab = args["export"] == "matlab";
// This part is export specific. For ALTA, the coefficients are just
......@@ -438,9 +438,9 @@ class compound_function: public nonlinear_function, public std::vector<nonlinear
// res += call_i(x);
for(int i=0; i<this->size(); ++i)
{
if(is_cpp || is_matlab || is_shader)
if(i != 0 && (is_cpp || is_matlab || is_shader))
{
out << "res += ";
out << "\tres += ";
}
this->at(i)->save_call(out, args);
......
#include "fitter.h"
#include <Eigen/Dense>
#include <ceres/ceres.h>
#include <string>
#include <iostream>
#include <fstream>
#include <limits>
#include <algorithm>
#include <cmath>
#include <QTime>
#include <core/common.h>
ALTA_DLL_EXPORT fitter* provide_fitter()
{
return new nonlinear_fitter_ceres();
}
class CeresFunctor : public ceres::CostFunction
{
public:
CeresFunctor(nonlinear_function* f, const vec xi) : _f(f), _xi(xi)
{
set_num_residuals(f->dimY());
mutable_parameter_block_sizes()->push_back(f->nbParameters());
}
virtual bool Evaluate(double const* const* x, double* y, double** dy) const
{
// Update the parameters vector
vec _p(_f->nbParameters());
for(int i=0; i<_f->nbParameters(); ++i) { _p[i] = x[0][i]; }
_f->setParameters(_p);
vec _di = vec(_f->dimY());
for(int i=0; i<_f->dimY(); ++i)
_di[i] = x[0][_f->dimX() + i];
// Should add the resulting vector completely
vec _y = _di - (*_f)(_xi);
for(int i=0; i<_f->dimY(); ++i)
y[i] = _y[i];
if(dy != NULL)
{
df(dy);
}
return true;
}
// The parameter of the function _f should be set prior to this function
// call. If not it will produce undesirable results.
void df(double ** fjac) const
{
// Get the jacobian of the function at position x_i for the current
// set of parameters (set prior to function call)
vec _jac = _f->parametersJacobian(_xi);
// For each output channel, update the subpart of the
// vector row
for(int i=0; i<_f->dimY(); ++i)
{
// Fill the columns of the matrix
for(int j=0; j<_f->nbParameters(); ++j)
{
fjac[0][i*_f->nbParameters() + j] = -_jac[i*_f->nbParameters() + j];
}
}
}
protected:
const vec _xi;
nonlinear_function* _f;
};
nonlinear_fitter_ceres::nonlinear_fitter_ceres()
{
}
nonlinear_fitter_ceres::~nonlinear_fitter_ceres()
{
}
bool nonlinear_fitter_ceres::fit_data(const data* d, function* fit, const arguments &args)
{
// I need to set the dimension of the resulting function to be equal
// to the dimension of my fitting problem
fit->setDimX(d->dimX()) ;
fit->setDimY(d->dimY()) ;
fit->setMin(d->min()) ;
fit->setMax(d->max()) ;
// Convert the function and bootstrap it with the data
if(dynamic_cast<nonlinear_function*>(fit) == NULL)
{
std::cerr << "<<ERROR>> the function is not a non-linear function" << std::endl;
return false;
}
nonlinear_function* nf = dynamic_cast<nonlinear_function*>(fit);
nf->bootstrap(d, args);
#ifndef DEBUG
std::cout << "<<DEBUG>> number of parameters: " << nf->nbParameters() << std::endl;
#endif
if(nf->nbParameters() == 0)
{
return true;
}
/* the following starting values provide a rough fit. */
vec p = nf->parameters();
// Create the problem
ceres::Problem problem;
for(int i=0; i<d->size(); ++i)
{
vec xi = d->get(i);
problem.AddResidualBlock(new CeresFunctor(nf, xi), NULL, &p[0]);
}
// Solver's options
ceres::Solver::Options options;
options.max_num_iterations = 25;
options.linear_solver_type = ceres::DENSE_QR;
options.minimizer_progress_to_stdout = true;
// Solves the NL problem
ceres::Solver::Summary summary;
ceres::Solve(options, &problem, &summary);
std::cout << summary.BriefReport() << std::endl;
std::cout << "<<INFO>> found parameters: " << p << std::endl;
nf->setParameters(p);
return true;
}
void nonlinear_fitter_ceres::set_parameters(const arguments& args)
{
}
#pragma once
// Include STL
#include <vector>
#include <string>
// Interface
#include <core/function.h>
#include <core/data.h>
#include <core/fitter.h>
#include <core/args.h>
#include <core/vertical_segment.h>
/*! \brief A fitter for non-linear BRDF models that uses Eigen's
* Levenberg-Marquardt solver.
* \ingroup plugins
*/
class nonlinear_fitter_ceres: public fitter
{
public: // methods
nonlinear_fitter_ceres() ;
virtual ~nonlinear_fitter_ceres() ;
// Fitting a data object
//
virtual bool fit_data(const data* d, function* fit, const arguments& args) ;
// Provide user parameters to the fitter
//
virtual void set_parameters(const arguments& args) ;
protected: // function
protected: // data
} ;
TARGET = nonlinear_fitter_ceres
TEMPLATE = lib
CONFIG *= qt \
plugin \
ceres \
eigen
DESTDIR = ../../build
INCLUDEPATH += ../..
HEADERS = fitter.h
SOURCES = fitter.cpp
LIBS += -L../../build \
-lcore
......@@ -58,12 +58,27 @@ void diffuse_function::load(std::istream &in)
void diffuse_function::save_call(std::ostream& out, const arguments& args) const
{
out << "#FUNC nonlinear_function_diffuse" << std::endl ;
for(int i=0; i<dimY(); ++i)
{
out << "kd " << _kd[i] << std::endl;
}
out << std::endl;
bool is_alta = !args.is_defined("export") || args["export"] == "alta";
if(is_alta)
{
out << "#FUNC nonlinear_function_diffuse" << std::endl ;
for(int i=0; i<dimY(); ++i)
{
out << "kd " << _kd[i] << std::endl;
}
out << std::endl;
}
else
{
out << "vec3(";
for(int i=0; i<dimY(); ++i)
{
out << _kd[i]; if(i < dimY()-1) { out << ", "; }
}
out << ")";
}
}
//! Number of parameters to this non-linear function
......
......@@ -378,7 +378,7 @@ void lafortune_function::load(std::istream& in)
{
in >> token >> _C[(n*_nY + i)*3 + 0];
in >> token >> _C[(n*_nY + i)*3 + 1];
in >> token >> _C[(n*_nY + i)*3 + 1];
in >> token >> _C[(n*_nY + i)*3 + 2];
in >> token >> _N[i];
}
......@@ -389,105 +389,75 @@ void lafortune_function::load(std::istream& in)
std::cout << "<<INFO>> N = " << _N << std::endl;
}
void lafortune_function::save(const std::string& filename) const
void lafortune_function::save_call(std::ostream& out, const arguments& args) const
{
std::ofstream file(filename.c_str(), std::ios_base::trunc);
file << "#DIM " << _nX << " " << _nY << std::endl ;
file << "#NB_LOBES " << _n << std::endl ;
for(int i=0; i<_nY; ++i)
{
file << _kd[i] << std::endl;
}
file << std::endl;
bool is_alta = !args.is_defined("export") || args["export"] == "alta";
for(int n=0; n<_n; ++n)
{
file << "#Lobe number " << n << std::endl;
file << "#Cx Cy Cz" << std::endl;
for(int i=0; i<_nY; ++i)
{
file << _C[(n*_nY + i)*3 + 0] << " " << _C[(n*_nY + i)*3 + 2] << " " << _C[(n*_nY + i)*3 + 1] << std::endl;
}
file << std::endl;
file << "#N" << std::endl;
for(int i=0; i<_nY; ++i)
{
file << _N[i] << std::endl;
}
file << std::endl;
}
}
//! \brief Output the function using a BRDF Explorer formating.
//! \todo Finish
void lafortune_function::save_brdfexplorer(const std::string& filename,
const arguments& args) const
{
std::ofstream file(filename.c_str(), std::ios_base::trunc);
file << "analytic" << std::endl;
file << std::endl;
// file << "::begin parameters" << std::endl;
// file << "::end parameters" << std::endl;
// file << std::endl;
file << std::endl;
file << "::begin shader" << std::endl;
type_definition(file, _nY) << " n;" << std::endl;
type_definition(file, _nY) << " Cx;" << std::endl;
type_definition(file, _nY) << " Cy;" << std::endl;
type_definition(file, _nY) << " Cz;" << std::endl;
type_definition(file, _nY) << " ";
type_affectation(file, std::string("kd"), _kd, _nY);
file << std::endl;
file << std::endl;
file << "vec3 BRDF( vec3 L, vec3 V, vec3 N, vec3 X, vec3 Y )" << std::endl;
file << "{" << std::endl;
if(_nY == 1)
if(is_alta)
{
file << " ";
type_definition(file, _nY) << " D = kd;" << std::endl << std::endl;
out << "#FUNC nonlinear_function_lafortune" << std::endl ;
out << "#NB_LOBES " << _n << std::endl ;
for(int n=0; n<_n; ++n)
{
file << " // Lobe number " << n+1 << std::endl;
file << " n = " << _N[n] << "; " << std::endl;
file << " Cx = " << _C[0*_n + n] << "; " << std::endl;
file << " Cy = " << _C[1*_n + n] << "; " << std::endl;
file << " Cz = " << _C[2*_n + n] << "; " << std::endl;
file << " D += pow(max(Cx * L.x * V.x + Cy * L.y * V.y + Cz * L.z * V.z, ";
type_definition(file, _nY) << "(0.0)), n);" << std::endl;
file << std::endl;
for(int i=0; i<_nY; ++i)
{
out << "Cx " << _C[(n*_nY + i)*3 + 0] << std::endl;
out << "Cx " << _C[(n*_nY + i)*3 + 1] << std::endl;
out << "Cz " << _C[(n*_nY + i)*3 + 2] << std::endl;
out << "N " << _N[n*_nY + i] << std::endl;
}
out << std::endl;
}
}
else
{
file << " ";
type_definition(file, _nY) << " D = kd;" << std::endl << std::endl;
for(int n=0; n<_n; ++n)
{
file << " // Lobe number " << n+1 << std::endl;
file << " "; type_affectation(file, std::string("n"), _N, _nY, n);
file << " "; type_affectation(file, std::string("Cx"), _C, _nY, n, 0, 3);
file << " "; type_affectation(file, std::string("Cy"), _C, _nY, n, 1, 3);
file << " "; type_affectation(file, std::string("Cz"), _C, _nY, n, 2, 3);
file << " D += pow(max(Cx * L.x * V.x + Cy * L.y * V.y + Cz * L.z * V.z, ";
type_definition(file, _nY) << "(0.0)), n);" << std::endl;
file << std::endl;
out << "lafortune(L, V, N, X, Y, vec3(";
for(int i=0; i<_nY; ++i)
{
out << _C[(n*_nY + i)*3 + 0];
if(i < _nY-1) { out << ", "; }
}
out << "), vec3(";
for(int i=0; i<_nY; ++i)
{
out << _C[(n*_nY + i)*3 + 1];
if(i < _nY-1) { out << ", "; }
}
out << "), vec3(";
for(int i=0; i<_nY; ++i)
{
out << _C[(n*_nY + i)*3 + 2];
if(i < _nY-1) { out << ", "; }
}
out << "), vec3(";
for(int i=0; i<_nY; ++i)
{
out << _N[n*_nY + i];
if(i < _nY-1) { out << ", "; }
}
// For multiple lobes, add a sum sign
out << "))";
if(n < _n-1) { out << " + "; }
}
}
file << " return vec3(D);" << std::endl;
// file << " if (normalized)" << std::endl;
// file << " D *= (2+n) / (2*PI);" << std::endl;
file << "}" << std::endl;
file << std::endl;
file << "::end shader" << std::endl;
file.close();
}
void lafortune_function::save_body(std::ostream& out, const arguments& args) const
{
out << "vec3 lafortune(vec3 L, vec3 V, vec3 N, vec3 X, vec3 Y, vec3 Cx, vec3 Cy, vec3 Cz, vec3 Nl)" << std::endl;
out << "{" << std::endl;
out << "\tvec3 ext_dot = Cx * dot(L,X)*dot(V,X + Cy * dot(L,Y)*dot(V,Y) + Cz * dot(L,N)*dot(V,N);" << std::endl;
out << "\treturn pow(max(ext_dot, vec3(0,0,0)), Nl);" << std::endl;
out << "}" << std::endl;
}
......@@ -28,9 +28,13 @@ class lafortune_function : public nonlinear_function
public: // methods
lafortune_function() : _n(1), _isotropic(false) { }
lafortune_function() : _n(1), _isotropic(false)
{
nonlinear_function::setParametrization(params::CARTESIAN);
nonlinear_function::setDimX(6);
}
// Overload the function operator
// Overload the function operator
virtual vec operator()(const vec& x) const ;
virtual vec value(const vec& x) const ;
virtual vec value(const vec& x, const vec& p);
......@@ -38,6 +42,10 @@ class lafortune_function : public nonlinear_function
//! \brief Load function specific files
virtual void load(std::istream& in) ;
// Save functions
void save_body(std::ostream& out, const arguments& args) const;
void save_call(std::ostream& out, const arguments& args) const;
//! \brief Boostrap the function by defining the diffuse term
virtual void bootstrap(const data* d, const arguments& args);
......@@ -90,15 +98,6 @@ class lafortune_function : public nonlinear_function
//! \brief Set the number of lobes to be used in the fit
void setNbLobes(int N);
protected: // methods
virtual void save(const std::string& filename) const;
//! \brief Output the function using a BRDF Explorer formating.
virtual void save_brdfexplorer(const std::string& filename,
const arguments& args) const;
private: // methods
//! \brief Provide the coefficient of the monochromatic lobe number
......
......@@ -10,6 +10,7 @@ SUBDIRS = \
rational_fitter_matlab \
# rational_fitter_dca \
nonlinear_levenberg_eigen \
nonlinear_fitter_ceres \
nonlinear_fresnel_schlick \
nonlinear_function_diffuse \
nonlinear_function_phong \
......
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