11d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// Ceres Solver - A fast non-linear least squares minimizer 21d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// Copyright 2013 Google Inc. All rights reserved. 31d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// http://code.google.com/p/ceres-solver/ 41d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// 51d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// Redistribution and use in source and binary forms, with or without 61d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// modification, are permitted provided that the following conditions are met: 71d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// 81d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// * Redistributions of source code must retain the above copyright notice, 91d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// this list of conditions and the following disclaimer. 101d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// * Redistributions in binary form must reproduce the above copyright notice, 111d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// this list of conditions and the following disclaimer in the documentation 121d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// and/or other materials provided with the distribution. 131d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// * Neither the name of Google Inc. nor the names of its contributors may be 141d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// used to endorse or promote products derived from this software without 151d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// specific prior written permission. 161d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// 171d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" 181d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE 191d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE 201d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE 211d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR 221d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF 231d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS 241d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN 251d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) 261d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE 271d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// POSSIBILITY OF SUCH DAMAGE. 281d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// 291d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// Author: mierle@gmail.com (Keir Mierle) 301d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// 311d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// An incomplete C API for Ceres. 321d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// 331d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// TODO(keir): Figure out why logging does not seem to work. 341d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 351d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling#include "ceres/c_api.h" 361d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 371d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling#include <vector> 381d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling#include <iostream> 391d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling#include <string> 401d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling#include "ceres/cost_function.h" 411d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling#include "ceres/loss_function.h" 421d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling#include "ceres/problem.h" 431d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling#include "ceres/solver.h" 441d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling#include "ceres/types.h" // for std 451d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling#include "glog/logging.h" 461d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 471d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingusing ceres::Problem; 481d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 491d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingvoid ceres_init() { 501d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // This is not ideal, but it's not clear what to do if there is no gflags and 511d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // no access to command line arguments. 52399f7d09e0c45af54b77b4ab9508d6f23759b927Scott Ettinger char message[] = "<unknown>"; 53399f7d09e0c45af54b77b4ab9508d6f23759b927Scott Ettinger google::InitGoogleLogging(message); 541d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling} 551d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 561d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingceres_problem_t* ceres_create_problem() { 571d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling return reinterpret_cast<ceres_problem_t*>(new Problem); 581d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling} 591d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 601d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingvoid ceres_free_problem(ceres_problem_t* problem) { 611d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling delete reinterpret_cast<Problem*>(problem); 621d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling} 631d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 641d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// This cost function wraps a C-level function pointer from the user, to bridge 651d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// between C and C++. 661d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingclass CallbackCostFunction : public ceres::CostFunction { 671d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling public: 681d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling CallbackCostFunction(ceres_cost_function_t cost_function, 691d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling void* user_data, 701d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling int num_residuals, 711d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling int num_parameter_blocks, 721d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling int* parameter_block_sizes) 731d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling : cost_function_(cost_function), 741d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling user_data_(user_data) { 751d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling set_num_residuals(num_residuals); 761d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling for (int i = 0; i < num_parameter_blocks; ++i) { 771d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling mutable_parameter_block_sizes()->push_back(parameter_block_sizes[i]); 781d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling } 791d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling } 801d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 811d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling virtual ~CallbackCostFunction() {} 821d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 831d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling virtual bool Evaluate(double const* const* parameters, 841d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling double* residuals, 851d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling double** jacobians) const { 861d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling return (*cost_function_)(user_data_, 871d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling const_cast<double**>(parameters), 881d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling residuals, 891d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling jacobians); 901d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling } 911d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 921d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling private: 931d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling ceres_cost_function_t cost_function_; 941d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling void* user_data_; 951d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling}; 961d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 971d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// This loss function wraps a C-level function pointer from the user, to bridge 981d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// between C and C++. 991d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingclass CallbackLossFunction : public ceres::LossFunction { 1001d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling public: 1011d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling explicit CallbackLossFunction(ceres_loss_function_t loss_function, 1021d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling void* user_data) 1031d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling : loss_function_(loss_function), user_data_(user_data) {} 1041d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling virtual void Evaluate(double sq_norm, double* rho) const { 1051d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling (*loss_function_)(user_data_, sq_norm, rho); 1061d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling } 1071d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 1081d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling private: 1091d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling ceres_loss_function_t loss_function_; 1101d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling void* user_data_; 1111d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling}; 1121d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 1131d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// Wrappers for the stock loss functions. 1141d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingvoid* ceres_create_huber_loss_function_data(double a) { 1151d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling return new ceres::HuberLoss(a); 1161d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling} 1171d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingvoid* ceres_create_softl1_loss_function_data(double a) { 1181d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling return new ceres::SoftLOneLoss(a); 1191d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling} 1201d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingvoid* ceres_create_cauchy_loss_function_data(double a) { 1211d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling return new ceres::CauchyLoss(a); 1221d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling} 1231d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingvoid* ceres_create_arctan_loss_function_data(double a) { 1241d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling return new ceres::ArctanLoss(a); 1251d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling} 1261d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingvoid* ceres_create_tolerant_loss_function_data(double a, double b) { 1271d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling return new ceres::TolerantLoss(a, b); 1281d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling} 1291d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 1301d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingvoid ceres_free_stock_loss_function_data(void* loss_function_data) { 1311d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling delete reinterpret_cast<ceres::LossFunction*>(loss_function_data); 1321d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling} 1331d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 1341d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingvoid ceres_stock_loss_function(void* user_data, 1351d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling double squared_norm, 1361d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling double out[3]) { 1371d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling reinterpret_cast<ceres::LossFunction*>(user_data) 1381d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling ->Evaluate(squared_norm, out); 1391d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling} 1401d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 1411d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingceres_residual_block_id_t* ceres_problem_add_residual_block( 1421d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling ceres_problem_t* problem, 1431d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling ceres_cost_function_t cost_function, 1441d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling void* cost_function_data, 1451d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling ceres_loss_function_t loss_function, 1461d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling void* loss_function_data, 1471d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling int num_residuals, 1481d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling int num_parameter_blocks, 1491d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling int* parameter_block_sizes, 1501d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling double** parameters) { 1511d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling Problem* ceres_problem = reinterpret_cast<Problem*>(problem); 1521d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 1531d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling ceres::CostFunction* callback_cost_function = 1541d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling new CallbackCostFunction(cost_function, 1551d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling cost_function_data, 1561d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling num_residuals, 1571d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling num_parameter_blocks, 1581d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling parameter_block_sizes); 1591d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 1601d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling ceres::LossFunction* callback_loss_function = NULL; 1611d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling if (loss_function != NULL) { 1621d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling callback_loss_function = new CallbackLossFunction(loss_function, 1631d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling loss_function_data); 1641d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling } 1651d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 1661d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling std::vector<double*> parameter_blocks(parameters, 1671d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling parameters + num_parameter_blocks); 1681d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling return reinterpret_cast<ceres_residual_block_id_t*>( 1691d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling ceres_problem->AddResidualBlock(callback_cost_function, 1701d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling callback_loss_function, 1711d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling parameter_blocks)); 1721d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling} 1731d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 1741d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingvoid ceres_solve(ceres_problem_t* c_problem) { 1751d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling Problem* problem = reinterpret_cast<Problem*>(c_problem); 176399f7d09e0c45af54b77b4ab9508d6f23759b927Scott Ettinger 1771d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // TODO(keir): Obviously, this way of setting options won't scale or last. 1781d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // Instead, figure out a way to specify some of the options without 1791d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // duplicating everything. 1801d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling ceres::Solver::Options options; 1811d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling options.max_num_iterations = 100; 1821d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling options.linear_solver_type = ceres::DENSE_QR; 1831d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling options.minimizer_progress_to_stdout = true; 1841d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 1851d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling ceres::Solver::Summary summary; 1861d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling ceres::Solve(options, problem, &summary); 1871d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling std::cout << summary.FullReport() << "\n"; 1881d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling} 189