Given parameter blocks $\left[x_{i_1}, \hdots , x_{i_k}\right]$, a
\texttt{CostFunction} is responsible for computing
a vector of residuals and if asked a vector of Jacobian matrices, i.e., given $\left[x_{i_1}, \hdots , x_{i_k}\right]$, compute the vector $f_i\left(x_{i_1},\hdots,x_{i_k}\right)$ and the matrices
\texttt{parameters} is an array of pointers to arrays containing the various parameter blocks. parameters has the same number of elements as parameter\_block\_sizes\_. Parameter blocks are in the same order as parameter\_block\_sizes\_.
\texttt{jacobians} is an array of size \texttt{parameter\_block\_sizes\_} containing pointers to storage for Jacobian matrices corresponding to each parameter block. The Jacobian matrices are in the same order as \texttt{parameter\_block\_sizes\_}. \texttt{jacobians[i]} is an array that contains \texttt{num\_residuals\_}$\times$\texttt{parameter\_block\_sizes\_[i]} elements. Each Jacobian matrix is stored in row-major order, i.e.,
If \texttt{jacobians} is \texttt{NULL}, then no derivatives are returned; this is the case when computing cost only. If \texttt{jacobians[i]} is \texttt{NULL}, then the Jacobian matrix corresponding to the $i^{\textrm{th}}$ parameter block must not be returned, this is the case when the a parameter block is marked constant.
\section{\texttt{SizedCostFunction}}
If the size of the parameter blocks and the size of the residual vector is known at compile time (this is the common case), Ceres provides \texttt{SizedCostFunction}, where these values can be specified as template parameters.
\begin{minted}{c++}
template<int kNumResiduals,
int N0 = 0, int N1 = 0, int N2 = 0, int N3 = 0, int N4 = 0, int N5 = 0>
In this case the user only needs to implement the \texttt{Evaluate} method.
\section{\texttt{AutoDiffCostFunction}}
But even defining the \texttt{SizedCostFunction} can be a tedious affair if complicated derivative computations are involved. To this end Ceres provides automatic differentiation.
To get an auto differentiated cost function, you must define a class with a
templated \texttt{operator()} (a functor) that computes the cost function in terms of
the template parameter \texttt{T}. The autodiff framework substitutes appropriate
\texttt{Jet} objects for T in order to compute the derivative when necessary, but
this is hidden, and you should write the function as if T were a scalar type
(e.g. a double-precision floating point number).
The function must write the computed value in the last argument (the only
non-\texttt{const} one) and return true to indicate success.
For example, consider a scalar error $e = k - x^\top y$, where both $x$ and $y$ are
two-dimensional vector parameters and $k$ is a constant. The form of this error, which is the
difference between a constant and an expression, is a common pattern in least
squares problems. For example, the value $x^\top y$ might be the model expectation
for a series of measurements, where there is an instance of the cost function
for each measurement $k$.
The actual cost added to the total problem is $e^2$, or $(k - x^\top y)^2$; however,
the squaring is implicitly done by the optimization framework.
To write an auto-differentiable cost function for the above model, first
define the object
\begin{minted}{c++}
class MyScalarCostFunction {
MyScalarCostFunction(double k): k_(k) {}
template <typename T>
bool operator()(const T* const x , const T* const y, T* e) const {
In the instantiation above, the template parameters following
\texttt{MyScalarCostFunction}, \texttt{<1, 2, 2>} describe the functor as computing a
1-dimensional output from two arguments, both 2-dimensional.
The framework can currently accommodate cost functions of up to 6 independent
variables, and there is no limit on the dimensionality of each of them.
\textbf{WARNING 1} Since the functor will get instantiated with different types for
\texttt{T}, you must convert from other numeric types to \texttt{T} before mixing
computations with other variables of type \texttt{T}. In the example above, this is
seen where instead of using \texttt{k\_} directly, \texttt{k\_} is wrapped with \texttt{T(k\_)}.
\textbf{WARNING 2} A common beginner's error when first using \texttt{AutoDiffCostFunction} is to get the sizing wrong. In particular, there is a tendency to
set the template parameters to (dimension of residual, number of parameters)
instead of passing a dimension parameter for {\em every parameter block}. In the
example above, that would be \texttt{<MyScalarCostFunction, 1, 2>}, which is missing
Here the convention is that the contribution of a term to the cost function is given by $\frac{1}{2}\rho(s)$, where $s =\|f_i\|^2$. Calling the method with a negative value of $s$ is an error and the implementations are not required to handle that case.
Most sane choices of $\rho$ satisfy:
\begin{align}
\rho(0) &= 0\\
\rho'(0) &= 1\\
\rho'(s) &< 1 \text{ in the outlier region}\\
\rho''(s) &< 0 \text{ in the outlier region}
\end{align}
so that they mimic the squared cost for small residuals.
\subsection{Scaling}
Given one robustifier $\rho(s)$
one can change the length scale at which robustification takes
place, by adding a scale factor $a > 0$ which gives us $\rho(s,a)= a^2\rho(s / a^2)$ and the first and second derivatives as $\rho'(s / a^2)$ and $(1/ a^2)\rho''(s / a^2)$ respectively.
\caption{Shape of the various common loss functions.}
\label{fig:loss}
\end{figure}
The reason for the appearance of squaring is that $a$ is in the units of the residual vector norm whereas $s$ is a squared norm. For applications it is more convenient to specify $a$ than
Here are some common loss functions implemented in Ceres. For simplicity we described their unscaled versions. Figure~\ref{fig:loss} illustrates their shape graphically.
where the terms involving the second derivatives of $f(x)$ have been ignored. Note that $H(x)$ is indefinite if $\rho''f(x)^\top f(x)+\frac{1}{2}\rho' < 0$. If this is not the case, then its possible to re-weight the residual and the Jacobian matrix such that the corresponding linear least squares problem for the robustified Gauss-Newton step.
In the case $2\rho''\left\|f(x)\right\|^2+\rho' \lesssim0$, we limit $\alpha\le1-\epsilon$ for some small $\epsilon$. For more details see Triggs et al~\cite{triggs-etal-1999}.
With this simple rescaling, one can use any Jacobian based non-linear least squares algorithm to robustifed non-linear least squares problems.
\texttt{GlobalSize} is the dimension of the ambient space in which the parameter block $x$ lives. \texttt{LocalSize} is the size of the tangent space that $\Delta x$ lives in. \texttt{Plus} implements $\boxplus(x,\Delta x)$ and $\texttt{ComputeJacobian}$ computes the Jacobian matrix
The \texttt{Problem} objects holds the robustified non-linear least squares problem~\eqref{eq:ceresproblem}. To create a least squares problem, use the \texttt{Problem::AddResidualBlock} and \texttt{Problem::AddParameterBlock} methods.
For example a problem containing 3 parameter blocks of sizes 3, 4 and 5
respectively and two residual blocks of size 2 and 6:
\texttt{AddResidualBlock} as the name implies, adds a residual block to the problem. It adds a cost function, an optional loss function, and connects the cost function to a set of parameter blocks.
The cost
function carries with it information about the sizes of the
parameter blocks it expects. The function checks that these match
the sizes of the parameter blocks listed in \texttt{parameter\_blocks}. The
program aborts if a mismatch is detected. \texttt{loss\_function} can be
\texttt{NULL}, in which case the cost of the term is just the squared norm
of the residuals.
The user has the option of explicitly adding the parameter blocks
using \texttt{AddParameterBlock}. This causes additional correctness
checking; however, \texttt{AddResidualBlock} implicitly adds the parameter
blocks if they are not present, so calling \texttt{AddParameterBlock}
Note that even though the Problem takes ownership of \texttt{cost\_function}
and \texttt{loss\_function}, it does not preclude the user from re-using
them in another residual block. The destructor takes care to call
delete on each \texttt{cost\_function} or \texttt{loss\_function} pointer only once,
regardless of how many residual blocks refer to them.
\texttt{AddParameterBlock} explicitly adds a parameter block to the \texttt{Problem}. Optionally it allows the user to associate a LocalParameterization object with the parameter block too. Repeated calls with the same arguments are ignored. Repeated
calls with the same double pointer but a different size results in undefined behaviour.
In fact you can set any number of parameter blocks to be constant, and Ceres is smart enough to figure out what part of the problem you have constructed depends on the parameter blocks that are free to change and only spends time solving it. So for example if you constructed a problem with a million parameter blocks and 2 million residual blocks, but then set all but one parameter blocks to be constant and say only 10 residual blocks depend on this one non-constant parameter block. Then the computational effort Ceres spends in solving this problem will be the same if you had defined a problem with one parameter block and 10 residual blocks.
\texttt{Problem} by default takes ownership of the
\texttt{cost\_function}, \texttt{loss\_function} and \\\texttt{local\_parameterization} pointers. These objects remain
live for the life of the \texttt{Problem} object. If the user wishes to
keep control over the destruction of these objects, then they can
do this by setting the corresponding enums in the \texttt{Options} struct. Even though \texttt{Problem} takes ownership of these pointers, it does not preclude the user from re-using them in another residual or parameter block. The destructor takes care to call