# Camera Motion Estimation from RGB-D-Inertial Scene Flow

Samuel Cerezo, Javier Civera
I3A, Universidad de Zaragoza
<EMAIL_ADDRESS>
###### Abstract

In this paper, we introduce a novel formulation for camera motion estimation
that integrates RGB-D images and inertial data through scene flow. Our goal is
to accurately estimate the camera motion in a rigid 3D environment, along with
the state of the inertial measurement unit (IMU). Our proposed method offers
the flexibility to operate as a multi-frame optimization or to marginalize
older data, thus effectively utilizing past measurements. To assess the
performance of our method, we conducted evaluations using both synthetic data
from the ICL-NUIM dataset and real data sequences from the OpenLORIS-Scene
dataset. Our results show that the fusion of these two sensors enhances the
accuracy of camera motion estimation when compared to using only visual data.

##  1 Introduction

Autonomous navigation plays a key role in enabling robots and various other
applications, including mixed reality and autonomous driving. For that,
precise motion estimates derived from onboard sensors are essential. And, in
this context, scene flow stands as one of the fundamental techniques for
motion estimation using RGB-D or range sensors [7, 9, 33, 30, 24]. More
specifically, scene flow refers to the estimation of the 3D motion field of
scene points obtained from two sensor readings [31]. Although optical and
scene flow have been used in numerous tasks over the years, such as motion
compensation [34], object tracking [21] and object learning [22], we focus in
this paper on the previously mentioned application of scene flow to odometry,
_i.e._ , the estimation of the camera motion. We are also motivated by
multisensor odometry and SLAM, which boost monocular-only approaches with
extra accuracy and robustness, highly relevant in safety-critical robotic
setups. Multisensor configurations have been widely explored in feature-based
odometry and SLAM, e.g., stereo cameras [20], visual-inertial [3] or LiDAR-
inertial [26], but are less explored in direct approaches that use the raw
sensor measurements without feature extraction.

On the one hand, RGB-D cameras provide a practical hardware alternative to
several challenges and limitations of visual odometry. Their availability at
low cost has facilitated many robotics and Augmented Reality (AR) applications
in the last decade. Today, RGB-D cameras stand out as one of the preferred
sensors for indoor applications in robotics and AR; and their future looks
promising either on their own or in combination with additional sensors. On
the other hand, most commercial mobile devices are equipped with Inertial
Measurement Units (IMU), which can provide large amount of information in
dynamic trajectories but exhibit large drift due to noises if not fused with
other information. This makes the visual-range-inertial fusion relevant, as
the three modalities offer complementary characteristics.

Our contribution in this paper is a RGB-D-inertial formulation for camera
motion estimation in rigid scenes. Up to our knowledge, this is the first time
that inertial data is fused together with color and depth measurements to
estimate camera motion based on optical flow. Specifically, we propose a
tightly coupled optimization by minimizing pre-integrated inertial residuals
and range constraints. As the inertial states are common between frames, we
formulate the problem as a multi-frame optimization, in which past frame’s
states can be estimated or marginalized out into prior residuals for the
inertial states. We evaluate our proposal in the synthetic ICL-NUIM dataset
and in the real OpenLORIS-Scene one. The effectiveness of our fusion is shown
by an error reduction of RGB-D-inertial estimation compared to RGB-D one.

##  2 Related Work

The first tracking system for ego-motion estimation which fuses vision and
inertial measurements was presented by Armesto _et al_. [1]. In this case, the
fusion is performed by considering a EKF and UKF (Extended and Unscented
Kalman Filters) with multi-rate sampling of measurements. The mentioned
sampling modality allows the system to work with the different rates of the
sensors. In 2013, Kerl _et al_. [10] proposed a fast and accurate method to
estimate the camera motion from RGB-D images. This approach estimates the
relative motion between two consecutive frames by minimizing the photometric
error. A motion prior is incorporated in the optimization, in order to guide
and stabilize motion estimation in the presence of dynamic objects. Nießner
_et al_. [17] developed an approach that improves the robustness of real-time
3D surface reconstruction by incorporating inertial sensing to the inter-frame
alignment. As a result, they could significantly reduce the number of
Iterative Closest Point (ICP) iterations required per frame. Modeling three-
dimensional scene motion as a twist field, Quiroga _et al_. [23] introduced a
method that encourages piecewise smooth solutions of rigid body motions. A
general formulation is given to solve local and global rigid motions by
jointly using intensity and depth data.

The first method to compute dense scene flow in real-time for RGB-D cameras
was introduced in 2015 by Jaimez _et al_. [8]. They proposed a variational
formulation where brightness and geometric consistency are imposed. Their
accuracy outperforms that of previous for RGB-D flow baselines, being able to
estimate non-rigid motions at 30Hz of frame rate. In the same year,
Leutenegger _et al_. [15] formulated a probabilistic cost function that
combines reprojection errors of landmarks with inertial terms, using stereo
and monocular cameras. On the other hand, a new dense method to compute the
odometry of a range sensor in real time is presented [7]. This method applies
the range flow constraint equation in order to obtain the velocity of the
sensor in a rigid environment. Experiments show that this approach
overperforms GICP which uses the same geometric input data, whereas it
achieves results similar to RDVO, which requires both geometric and
photometric data.

The first tightly-coupled dense RGB-D-inertial SLAM system was proposed in
2017 by Laidlow _et al_. [13]. This system jointly optimises the camera pose,
velocity, IMU biases and gravity direction while building up a globally
consistent, fully dense surfel-based 3D reconstruction of the environment. In
2019, Shan _et al_. introduced VINS-RGBD [27]. The authors integrate a mapping
system based on depth data and octree filtering to achieve real-time mapping.
However the proposed system is applied only in ground robots. A RGB-D scene
flow estimation method with global nonrigid and local rigid motion assumption
is proposed by Li _et al_. in [16]. 3D motion is estimated based on the global
non-rigid and local rigid assumption and spatial-temporal correlation of RGBD
information. With this assumption, the interaction of motion from different
parts in the same segmented region is avoided, especially the non-rigid
object, e.g., a human body.

The flow formulation has been adapted to novel sensor modalities, e.g., event
cameras [25], or to include additional information, such as robot dynamics in
the work of Lee _et al_. [14]. Similarly to ours, the motivation in this last
case is improving the robustness and accuracy of the camera motion estimation.
Zhai _et al_. [32] compiled in a survey recent advances on optical and scene
flow. Up to our knowledge, inertial sensing has never been integrated in flow
formulations. Our work contributes to the literature presenting the first
camera motion estimation from RGB-D-inertial scene flow, demonstrating its
effectiveness in simulated and real public datasets.

##  3 Notation

Throughout this article, bold lower-case letters ($\mathbf{x}$) represent
vectors and bold upper-case letters ($\boldsymbol{\Sigma}$) matrices. Scalars
will be represented by light lower-case letters (${\alpha}$), scalar functions
and images by light upper-case letters ($J$). Camera poses are represented as
$\mathbf{T}_{WB}=[\mathbf{R}_{WB},\hskip 1.42262pt^{W}\mathbf{p}]\in{SE}(3)$
and transform points from frame $B$ to world coordinate system $W$.

##  4 IMU Model and Motion Integration

###  4.1 Inertial preintegration

An IMU consists typically of an accelerometer and a three-axis gyroscope, and
measures the angular velocity ${}^{B}\boldsymbol{\omega}$ and linear
acceleration ${}^{B}\mathbf{a}$ of the sensor in the body reference frame $B$.
We will denote the IMU measurement at time $k$ as
${}^{B}\tilde{\boldsymbol{\omega}}_{k}$ and ${}^{B}\tilde{\mathbf{a}}_{k}$.
IMU measurements are affected by additive white noise $\boldsymbol{\eta}^{g}$,
$\boldsymbol{\eta}^{a}\in\mathbb{R}^{3}$ and two slowly varying gyroscope and
accelerometer bias $\mathbf{b}^{g}$ and $\mathbf{b}^{a}\in\mathbb{R}^{3}$
respectively. Finally, the acceleration measurement is affected by gravity
${}^{W}\mathbf{g}$. This model is formulated by Eq. (1) and (2).

$^{B}\tilde{\boldsymbol{\omega}}_{k}=\hskip
1.42262pt^{B}\boldsymbol{\omega}_{k}+\mathbf{b}^{g}_{k}+\boldsymbol{\eta}^{g}_{k}$
(1)
$^{B}\tilde{\mathbf{a}}_{k}=\mathbf{R}_{WB}^{\top}\left({}^{W}\mathbf{a}_{k}-\hskip
1.42262pt^{W}\mathbf{g}\right)+\mathbf{b}^{a}_{k}+\boldsymbol{\eta}^{a}_{k}$
(2)

We use pre-integrated inertial residuals as proposed by Forster _et al._ [4].
We compute an inital guess for $\mathbf{b}^{g}$ as the difference of an
estimate ${}^{B}\boldsymbol{\omega}_{k}$ of the angular velocity between two
consecutive frames, and the direct measurement of the gyroscope
${}^{B}\tilde{\boldsymbol{\omega}}_{k}$ that includes the bias, _i.e._ ,
$\hat{\mathbf{b}}^{g}_{k}=\hskip
1.42262pt^{B}\tilde{\boldsymbol{\omega}}_{k}-\hskip
1.42262pt^{B}\boldsymbol{\omega}_{k}$, where the angular velocity estimate
${}^{B}\boldsymbol{\omega}_{k}$ can be computed by relative motion estimation
between point clouds, divided by the time increment between them. We set the
initial seed for $\mathbf{b}^{a}_{k}$ to zero.

Following [4], we can use the relative motion increment
$\Delta\mathbf{v}_{ij}\doteq\hskip
1.42262pt\mathbf{R}_{i}^{\top}(\mathbf{v}_{j}-\mathbf{v}_{i}-\mathbf{g}\Delta
t_{ij})$ in order to obtain a first gravity vector estimation

$\hat{\mathbf{g}}=\frac{\mathbf{v}_{j}-\mathbf{v}_{i}}{\Delta
t_{ij}}-\frac{\mathbf{R}_{i}\Delta\mathbf{v}_{ij}}{\Delta t_{ij}}$ (3)

###  4.2 Noise propagation

The covariance matrix of the raw IMU measurements noise
$\boldsymbol{{\Sigma}_{\eta}}\in\mathbb{S}^{6}_{+}$111By
$\mathbb{S}^{n}_{+}=\\{\mathbf{\Sigma}\in\mathbb{R}^{n\times n}\ |\
\mathbf{\Sigma}=\mathbf{\Sigma}^{\top},\mathbf{\Sigma}\succeq 0\\}$ we denote
the set of $n\times n$ symmetric positive semidefinite matrices. is composed
by sub-matrices
$\boldsymbol{{\Sigma}_{\omega}},\boldsymbol{{\Sigma}_{\textbf{a}}}\in\mathbb{S}^{3}_{+}$

$\boldsymbol{{\Sigma}_{\eta}}=\begin{bmatrix}\begin{array}[]{c:c}\boldsymbol{{\Sigma}_{\omega}}&\textbf{0}_{3\mathrm{x}3}\\\
\hdashline\textbf{0}_{3\mathrm{x}3}&\boldsymbol{{\Sigma}_{\textbf{a}}}\end{array}\end{bmatrix}$
(4)

where $\textbf{0}_{3\mathrm{x}3}$ is a $3\mathrm{x}3$ matrix for which all its
elements are equal to zero.

Following the computation of the preintegrated noise covariance in [4], we
consider the matrix
$\boldsymbol{\eta}_{ik}^{\Delta}\doteq\left[\delta\boldsymbol{\phi}_{ik}^{\top},\delta\mathbf{v}_{ik}^{\top},\delta\mathbf{p}_{ik}^{\top}\right]^{\top}$
and defining the IMU measurement noise
$\boldsymbol{\eta}_{k}^{d}\doteq\left[\boldsymbol{\eta}_{k}^{gd},\boldsymbol{\eta}_{k}^{ad}\right]^{\top}$,
the noise is propagated as

${\boldsymbol{\Sigma}}_{ij}=\mathbf{A}_{j-1}{\boldsymbol{\Sigma}}_{ij-1}\mathbf{A}_{j-1}^{\top}+\mathbf{B}_{j-1}\boldsymbol{{\Sigma}_{\eta}}\mathbf{B}_{j-1}^{\top}$
(5)

with initial conditions ${\boldsymbol{\Sigma}}_{ii}=\textbf{0}_{9\times 9}$
and $\mathbf{A}_{j-1}\in\mathbb{R}^{9\times 9}$ ,
$\mathbf{B}_{j-1}\in\mathbb{R}^{9\times 6}$ defined as

$\mathbf{A}_{j-1}=\begin{bmatrix}\Delta\tilde{\mathbf{R}}_{j-1j}^{\top}&\textbf{0}_{3\mathrm{x}3}&\textbf{0}_{3\mathrm{x}3}\\\
-\Delta\tilde{\mathbf{R}}_{ij-1}(\tilde{\mathbf{a}}_{j-1}-\mathbf{b}_{i}^{a})^{\wedge}\Delta
t&\textbf{I}_{3}&\textbf{0}_{3\mathrm{x}3}\\\
-\frac{1}{2}\Delta\tilde{\mathbf{R}}_{ij-1}(\tilde{\mathbf{a}}_{j-1}-\mathbf{b}_{i}^{a})^{\wedge}\Delta
t^{2}&\textbf{I}_{3}\Delta t&\textbf{I}_{3}\end{bmatrix}$
$\mathbf{B}_{j-1}=\begin{bmatrix}\mathbf{J}_{r}^{j-1}\Delta
t&\textbf{0}_{3\mathrm{x}3}\\\
\textbf{0}_{3\mathrm{x}3}&\Delta\tilde{\mathbf{R}}_{ij-1}\Delta t\\\
\textbf{0}_{3\mathrm{x}3}&\frac{1}{2}\Delta\tilde{\mathbf{R}}_{ij-1}\Delta
t^{2}\end{bmatrix}$

where $\textbf{I}_{3}$ stands for the identity matrix of size $3$.

For the computation of the above matrices, the preintegrated expressions in
[4] were used. The matrix ${\boldsymbol{\Sigma}}_{ij}\in\mathbb{S}^{9}_{+}$ is
composed by nine sub-matrices of dimension 3 by 3 each

$\boldsymbol{{\Sigma}}_{ij}=\begin{bmatrix}\begin{array}[]{c:c:c}{\boldsymbol{\Sigma}}_{\Delta\boldsymbol{\phi}_{ij}}&{\boldsymbol{\Sigma}}_{\Delta\boldsymbol{\phi}_{ij}\Delta\mathbf{v}_{i}}&{\boldsymbol{\Sigma}}_{\Delta\boldsymbol{\phi}_{ij}\Delta\mathbf{p}_{ij}}\\\
\hdashline{\boldsymbol{\Sigma}}_{\Delta\mathbf{v}_{i}\Delta\boldsymbol{\phi}_{ij}}&{\boldsymbol{\Sigma}}_{\Delta\mathbf{v}_{i}}&{\boldsymbol{\Sigma}}_{\Delta\mathbf{v}_{i}\Delta\mathbf{p}_{ij}}\\\
\hdashline{\boldsymbol{\Sigma}}_{\Delta\mathbf{p}_{ij}\Delta\boldsymbol{\phi}_{ij}}&{\boldsymbol{\Sigma}}_{\Delta\mathbf{p}_{ij}\Delta\mathbf{v}_{i}}&{\boldsymbol{\Sigma}}_{\Delta\mathbf{p}_{ij}}\end{array}\end{bmatrix}$
(6)

The above matrix will be important in the calculation of inertial residuals in
the next section.

###  4.3 Gravity vector representation

As the gravity modulus is known, a reasonable representation is by its
directional vector. Unit-norm direction vectors belong to the ${S}^{2}$
manifold, which has only two degrees of freedom. As ${S}^{2}$ does not form a
Lie group, we follow the parametrization proposed in [6]. The explicit
expressions are detailed in [19].

###  4.4 Optical flow and velocity constraint

We adopt the standard assumption that the local intensity image patterns are
approximately constant, at least in the short period of time between two
frames of a video [2]. This constraints the motion in the image as in the
following

$0=\frac{\partial I}{\partial t}+\frac{\partial I}{\partial
u}{\dot{u}}+\frac{\partial I}{\partial v}{\dot{v}}$ (7)

where $(\dot{u},\dot{v})$ is the optical flow in image units
($\text{pixels}/s$). Under the common assumption of a rigid scene, we
formulate the point velocities in terms of the camera motion. Let
$Z:\Omega\rightarrow\mathbb{R}$ be a depth image provided by a 3D range camera
where $\Omega$ is the image domain. Following the work by Spies et al. [29],
the range flow constraint is as follows

${\dot{Z}}=({\dot{w}})=\frac{\partial Z}{\partial t}+\frac{\partial
Z}{\partial u}{\dot{u}}+\frac{\partial Z}{\partial v}{\dot{v}}$ (8)

This equation reflects that the total derivative of the depth can be
calculated from the optical flow and the partial derivatives of $Z$. Following
[7] and using the pin-hole model, we obtain the range flow constraint in Eq.
(4.4). Here $f_{x},f_{y}$ are the focal length values, expressed in pixels,
while $\dot{x},\dot{y},\dot{z}$ are in camera coordinates.

$\displaystyle-\frac{\partial Z}{\partial t}$ $\displaystyle=$
$\displaystyle\left(1+\frac{xf_{x}}{z^{2}}\frac{\partial Z}{\partial
u}+\frac{yf_{y}}{z^{2}}\frac{\partial Z}{\partial
v}\right)\left(v_{z}+y\omega_{x}-x\omega_{y}\right)$
$\displaystyle+\frac{f_{x}}{z}\frac{\partial Z}{\partial
u}(-v_{x}+y\omega_{z}-z\omega_{y})$
$\displaystyle+\frac{f_{y}}{z}\frac{\partial Z}{\partial
v}(-v_{y}-x\omega_{z}+z\omega_{x})$

The above constraint for the camera velocity will be used for our visual
residuals, detailed in next section.

##  5 Camera Motion from RGB-D-I Flow

This section presents our approach to integrating inertial measurements with
RGB-D scene flow to estimate camera motion. We begin by defining the state,
which varies according to different operating modes, primarily depending on
the number of frames considered. We then proceed to formulate the cost
function to be optimized. Finally, this section concludes with the
marginalization process, through which we retain the information of removed
states.

###  5.1 State definition

Our goal is to track the state $\mathbf{x}$ of a sensing device equipped with
an IMU and a RGB-D camera. This state consists basically of the device
velocities, IMU biases and gravity vector at different moments of time. We
assume that the IMU is synchronized with the camera, as it is shown in Fig. 1.

Figure 1: Illustration of the temporal notation for RGB-D images, IMU
measurements and marginalization and optimization windows.

In a general case, the system performs the optimization over N frames, what we
called a N-sliding window. As the sensing device moves through the trajectory,
we marginalize out the old states for the optimization to remain compact. The
state at step $k+N-1$ is defined as

$\mathbf{x}=\left[{\mathbf{v}}_{k}^{\top},{\boldsymbol{\omega}}_{k}^{\top},\cdots,{\mathbf{v}}_{k+N-1}^{\top},{\boldsymbol{\omega}}_{k+N-1}^{\top},\mathbf{g}^{\top},{\mathbf{b}^{g}}^{\top},{\mathbf{b}^{a}}^{\top}\right]^{\top}$
(10)

where $\\{\mathbf{v}_{l}\\}^{k+N-1}_{l=k}\in\mathbb{R}^{3}$ are the linear
velocities in each frame,
$\\{\boldsymbol{\omega}_{l}\\}^{k+N-1}_{l=k}\in\mathbb{R}^{3}$ are the angular
velocities in each frame, $\mathbf{g}\in\mathrm{S}^{2}$ contain the two
degrees of freedom of the gravity direction and
$\mathbf{b}^{g},\mathbf{b}^{a}\in\mathbb{R}^{3}$ are the gyroscope and
accelerometer bias, respectively.

As mentioned, we include the gravity direction as part of the state. Under the
traditional absolute formulations for visual-inertial state estimation, this
variable is removed from the state by aligning the global reference frame to
the gravity direction during system initialization. But, as a consequence, the
rest of the variables in the state are tied to a gravity-aligned absolute
frame. In turn, by including the gravity vector in the local camera frames we
remove this dependence and all states are relative. As an additional benefit,
it becomes possible to explicitly re-estimate the gravity direction during
normal system operation and thus avoid coupling gravity and absolute
orientation errors. In order to improve the observability of the state (in
particular of the accelerometer bias), we assume a known gravity magnitude
(981 $cm/s^{2}$) and only optimize the gravity direction.

###  5.2 Cost function

We formulate an optimization problem over the state $\mathbf{x}$ for which the
camera velocity consistency is imposed as well as those terms corresponding to
the pre-integration of the IMU readings. The joint optimization problem will
consist on minimizing a cost function $J(\mathbf{x})$ which is the summation
of terms associated to the inertial measurements $J_{i}$ as well as to the
camera measurements $J_{c}$. Our state estimate $\hat{\mathbf{x}}$ will be the
one that minimizes the cost function $J(\mathbf{x})$.

$\displaystyle\hat{\mathbf{x}}=\operatorname*{arg\hskip 1.42262ptmin\hskip
2.84526pt}_{\mathbf{x}}J(\mathbf{x})=\operatorname*{arg\hskip
1.42262ptmin\hskip
2.84526pt}_{\mathbf{x}}\left(J_{c}(\mathbf{x})+J_{i}(\mathbf{x})\right)$ (11)

We begin by developing the term $J_{c}$. For a pair of consecutive frames $i$
and $j$, the velocity constraint in Eq. (4.4) results in the linear constraint

$\mathbf{r}_{c}=\mathbf{WAx}-\mathbf{WB}$ (12)

where $\mathbf{A}$ contains the weights of the coefficients that multiply the
state vector $\mathbf{x}$ in the velocity constraint, and the matrix
$\mathbf{B}$ contains the temporal derivatives of the per-pixel depths
(inverted in sign). The linearization that is applied to derive the range flow
constraint in Eq. (4.4) assumes differentiability of the depth images and
small scene displacement. Therefore, we implement an adaptive mask on the
image, in order to discard those pixels belonging to edges and prone to have
high depth derivatives. This mask is represented in a diagonal matrix
$\mathbf{W}$, which also has the weights associated with the uncertainty of
each equation. For details on these aspects, the reader is referred to [7].
Using the residual in Eq. (12), the visual cost is expressed as

$J_{c}=\mathbf{r}_{c}^{\top}{\boldsymbol{\Sigma}}_{c}^{-1}\mathbf{r}_{c}$ (13)

Having developed the first term of the cost function, we now turn our
attention to the development of $J_{i}$, which will be composed by several
terms. Firstly, from Eq. (1) we can derive the residual associated to the
angular velocity estimate, as follows

$\mathbf{r}_{\omega}=\hskip
1.42262pt^{B}\boldsymbol{\omega}-\left({}^{B}\tilde{\boldsymbol{\omega}}-\mathbf{b}^{g}\right)$
(14)

Assuming constant biases, as we compute flow for a small number of frames, the
residual for the preintegrated linear velocity term is defined in [4] as

$\mathbf{r}_{\Delta\mathbf{v}_{i}}=\hskip
1.42262pt\mathbf{R}_{i}^{\top}(\mathbf{v}_{j}-\mathbf{v}_{i}-\mathbf{g}\Delta
t)-\Delta\tilde{\mathbf{v}}_{ij}$ (15)

The residual for biases are made by penalising changes as the sliding window
moves _i.e._ , $\mathbf{r}_{{ba}}=\mathbf{b}_{0}^{a}-\mathbf{b}^{a}$ and
$\mathbf{r}_{{bg}}=\mathbf{b}_{0}^{g}-\mathbf{b}^{g}$, where
$\mathbf{b}_{0}^{a}$ and $\mathbf{b}_{0}^{g}$ are the initial estimates.

Using the submatrices
${\boldsymbol{\Sigma}}_{\Delta\mathbf{v}_{i}},{\boldsymbol{\Sigma}}_{\boldsymbol{\omega}},{\boldsymbol{\Sigma}}_{\mathbf{a}}\in\mathbb{S}^{3}_{+}$
from Eq. (6), we can now define the term $J_{i}$ as follows

$\displaystyle J_{i}$ $\displaystyle=$
$\displaystyle\mathbf{r}_{\Delta\mathbf{v}_{i}}^{\top}{\boldsymbol{\Sigma}}_{\Delta\mathbf{v}_{i}}^{-1}\mathbf{r}_{\Delta\mathbf{v}_{i}}+\mathbf{r}_{\omega}^{\top}{\boldsymbol{\Sigma}}_{\boldsymbol{\omega}}^{-1}\mathbf{r}_{\omega}$
(16)
$\displaystyle+\mathbf{r}_{{bg}}^{\top}{\boldsymbol{\Sigma}}_{\boldsymbol{\omega}}^{-1}\mathbf{r}_{{bg}}+\mathbf{r}_{{ba}}^{\top}{\boldsymbol{\Sigma}}_{\mathbf{a}}^{-1}\mathbf{r}_{{ba}}$

Up to this point we have defined the cost function, by means of $J_{c}$ and
$J_{i}$. In the general case $J(\mathbf{x})$ will be made up depending on the
number of frames in each case. These cases will be detailed below.

###  5.3 Operating Modes

We have mentioned that the cost function $J(\mathbf{x})$ will be formed as a
function of the number of frames ($N$) in the sliding window while the camera
is moving along the trajectory. Fig. 2 illustrates this situation graphically.

Figure 2: Optimization using a sliding window over N frames while the camera
is moving over the trajectory.

Depending on the value of $N$, the system will have different operating modes
which will be detailed in this section.
We start with the factor graph illustrated in Fig. 3(a), where only frames $i$
and $j$ are available. Here, $N=2$ and any variable has been marginalized yet.
Since we only have two frames, the cost function will be
$J=J_{i}^{ij}+J_{c}^{ij}$. The first term is associated to the inertial
measurements and the second one to the camera measurements. Super-indices
${}^{ij}\cdot$ denote that the corresponding term is built up by frames $i$
and $j$. The state $\mathbf{x}$ contains in this case the velocities
${\mathbf{v}}_{i},{\boldsymbol{\omega}}_{i},{\mathbf{v}}_{j}$ and
${\boldsymbol{\omega}}_{j}$ as well as the gravity ${}^{W}\mathbf{g}$ and
biases $\mathbf{b}^{g}$ and $\mathbf{b}^{a}$.

Figure 3: Factor graph representation using different modes of operation. Blue
and green shapes contains the variables to be estimated. (a) Taking two
frames, only one visual and one inertial residual are used. (b) Here we take
three frames, so there are two visual residual. Only one inertial residual is
used. (c) Is the same situation as before but in this case two inertial
residual are used. The difference between (b) and (c) is the inertial
constraint imposed by the last aggregate frame.

When a new frame ($k-$frame) is added to the window, $N=3$ and the associated
graph is shown in Fig. 3(b). In this case the state remains the same as before
but the cost changes to $J(\mathbf{x})=J_{i}^{ij}+J_{c}^{ij}+J_{c}^{jk}$,
_i.e._ , only visual information is added to the cost function.
The last case is when we add a new inertial term (between $j$ and $k$ frames),
therefore $J(\mathbf{x})=J_{i}^{ij}+J_{c}^{ij}+J_{c}^{jk}+J_{i}^{jk}$. This
situation is represented in Fig. 3(c). Here the state changes and the
velocities ${\mathbf{v}}_{k}$ and ${\boldsymbol{\omega}}_{k}$ are added.

In the general case, the state $\mathbf{x}\in\mathbb{R}^{6N+8}$ is defined as

$\mathbf{x}=\left[{\mathbf{v}}_{i}^{\top},{\boldsymbol{\omega}}_{i}^{\top},\dots,{\mathbf{v}}_{i+N-1}^{\top},{\boldsymbol{\omega}}_{i+N-1}^{\top},{\mathbf{g}}^{\top},{\mathbf{b}^{g}}^{\top},{\mathbf{b}^{a}}^{\top}\right]^{\top}$
(17)

and the cost function $J(\mathbf{x})$ can be expressed compactly as follows

$J(\mathbf{x})=\sum_{p=i}^{i+N-1}\left(\mathbf{r}_{c_{p}}^{\top}{\boldsymbol{\Sigma}}_{c_{p}}^{-1}\mathbf{r}_{c_{p}}+\mathbf{r}_{\Delta\mathbf{v}_{p}}^{\top}{\boldsymbol{\Sigma}}_{\Delta\mathbf{v}_{p}}^{-1}\mathbf{r}_{\Delta\mathbf{v}_{p}}\right)+\mathbf{r}_{bg}^{\top}{\boldsymbol{\Sigma}}_{\boldsymbol{\omega}}^{-1}\mathbf{r}_{bg}+\mathbf{r}_{ba}^{\top}{\boldsymbol{\Sigma}}_{\mathbf{a}}^{-1}\mathbf{r}_{ba}+\sum_{l=i}^{i+N}\mathbf{r}_{\omega_{l}}^{\top}{\boldsymbol{\Sigma}}_{\boldsymbol{\omega}}^{-1}\mathbf{r}_{\omega_{l}}$
(18)

###  5.4 Marginalization

We mentioned that, as the sliding window moves, information from previous
states is marginalized out. Let consider the case in Fig. 4, in which we want
to perform 3-frames-sliding-window optimization.

Figure 4: Factor graph using sliding window with containing 3 frames. Blue and
green shapes contains the variables to be estimated. (a) When a new frame
comes, both visual and inertial residual is added and the marginalization is
done. (b) After marginalization, a new prior residual is added on the cost
function.

Fig. 4(a) illustrates the graph after we have optimized the states at time
steps $i$, $j$ and $k$, and a new frame ($l-$frame) arrives. The state
estimate at this point, ${\mathbf{x}}_{M}$, can be decomposed into two parts:
the variables that we want to marginalize
$\mathbf{x}_{\alpha}=[\mathbf{v}_{i}^{\top},\boldsymbol{\omega}_{i}^{\top}]^{\top}\in\mathbb{R}^{6}$
and the starting point for the next step
${\mathbf{x}}_{\beta}=\left[{\mathbf{v}}_{j}^{\top},{\boldsymbol{\omega}}_{j}^{\top},{\mathbf{v}}_{k}^{\top},{\boldsymbol{\omega}}_{k}^{\top},{\mathbf{v}}_{l}^{\top},{\boldsymbol{\omega}}_{l}^{\top},{\mathbf{g}}^{\top},{\mathbf{b}^{g}}^{\top},{\mathbf{b}^{a}}^{\top}\right]^{\top}$.

${\mathbf{x}}_{M}=\begin{bmatrix}\begin{array}[]{c:c}\mathbf{x}_{\alpha}^{\top}&\mathbf{x}_{\beta}^{\top}\end{array}\end{bmatrix}^{\top}$
(19)

Consider Eq. (20), where we have the cost function before marginalization
$J(\mathbf{x})$, and a new term $J_{p}$, which stands for the marginalization
priors and accounts for the information associated to the marginalized
variables

$J^{*}(\mathbf{x})\doteq J(\mathbf{x})+J_{p}$ (20)

Using second-order Taylor approximation, the cost $J(\mathbf{x})$ can be
expressed as follows

$J(\mathbf{x})\approx J(\mathbf{x}_{0})+\cancel{\nabla
J(\mathbf{x}_{0})\mathbf{r}}+{\frac{1}{2}\mathbf{r}^{\top}\mathbf{H}(\mathbf{x}_{0})\mathbf{r}}$
(21)

The above approximation is calculated around a state $\mathbf{x}_{0}$, where a
minimum is achieved, _i.e._ $\nabla J(\mathbf{x}_{M})\mathbf{r}=0$
($\mathbf{r}=\mathbf{x}-\mathbf{x}_{0}$). The Hessian $\mathbf{H}$ contains
the second derivatives of the cost function with respect to the state
variables, therefore it encodes how every state variable affects the others.

We denote as $\alpha$ the block of variables we would like to marginalize, and
$\beta$ the block of variables we would like to keep. When marginalizing a set
$\alpha$ of variables, we gather all factors dependent on them as well as the
connected variables $\beta$. This is done by means of the Schur Complement as
follows

$\mathbf{H}^{*}=\mathbf{H}_{\beta\beta}-\mathbf{H}_{\alpha\beta}^{\top}\mathbf{H}_{\alpha\alpha}^{-1}\mathbf{H}_{\alpha\beta}$
(22)

Fig. 5 graphically illustrates how the $\alpha-$block (variables
$\mathbf{v}_{i}$ and $\boldsymbol{\omega}_{i}$) is removed from $\mathbf{H}$
but the information is preserved in the new matrix $\mathbf{H}^{*}$.

Figure 5: Marginalization example. We start with a Hessian matrix $\mathbf{H}$
after optimization with $N=4$. We want to marginalize $\mathbf{v}_{i}$ and
$\boldsymbol{\omega}_{i}$. The marginalized Hessian matrix $\mathbf{H}^{*}$
corresponds to the Schur complement of $\mathbf{H}_{\alpha\alpha}$. This
calculation transfers the information constraints of the variable being
eliminated to its adjacent nodes, adding shared information between these
variables (green cells).

Considering ${\mathbf{x}}_{\beta}$ obtained from Eq. (19) and the new state
$\mathbf{x}^{*}$, the term $\mathbf{r}_{p}$ can be expressed as
$\mathbf{r}_{p}=\mathbf{x}^{*}-{\mathbf{x}}_{\beta}$, where

$\mathbf{x}^{*}=\left[{\mathbf{v}}_{j}^{\top},{\boldsymbol{\omega}}_{j}^{\top},{\mathbf{v}}_{k}^{\top},{\boldsymbol{\omega}}_{k}^{\top},{\mathbf{v}}_{l}^{\top},{\boldsymbol{\omega}}_{l}^{\top},{\mathbf{g}}^{\top},{\mathbf{b}^{g}}^{\top},{\mathbf{b}^{a}}^{\top}\right]^{\top}$

The new term $J_{p}$ can now be defined as follows

$J_{p}=\mathbf{r}_{p}^{\top}\hskip 0.85358pt{\mathbf{H}^{*}}\hskip
0.85358pt\mathbf{r}_{p}$ (23)

Finally, the result is obtained by minimizing the cost function
$J^{*}(\mathbf{x})$ expressed in Eq. (20).

| 2-frames | 3-frames | 4-frames | 5-frames
---|---|---|---|---
|  | RGB-D
---
(DIFODO∗)
| RGB-D-I
---
(ours)
| RGB-D
---
(∗∗)
| RGB-D-I
---
(ours)
| RGB-D-I
---
(ours+M)
| RGB-D
---
(∗∗)
| RGB-D-I
---
(ours)
| RGB-D
---
(∗∗)
| RGB-D-I
---
(ours)
| RGB-D-I
---
(ours+M)
| $\mathrm{RMSE}_{{\mathbf{v}}}$
---
[$cm/s$]
| 8.020
---
±4.989
| 7.897
---
±4.699
| 0.565
---
±0.508
| 0.524
---
±0.442
| 0.523
---
±0.250
| 0.575
---
±0.514
| 0.533
---
±0.447
| 0.582
---
±0.520
| 0.538
---
±0.452
| 0.535
---
±0.245
| $\mathrm{RMSE}_{{\boldsymbol{\omega}}}$
---
[$rad/s$]
| 0.168
---
±0.079
| 0.037
---
±0.022
| 0.00113
---
±0.00081
| 0.00107
---
±0.00072
| 0.00091
---
±0.00052
| 0.00114
---
±0.00083
| 0.00108
---
±0.00074
| 0.00116
---
±0.00085
| 0.00109
---
±0.00075
| 0.00094
---
±0.00054
| $\mathrm{RMSE}_{\mathbf{b}^{a}}$
---
[$cm/s^{2}$]
- |  | 0.177
---
±0.026
- |  | 0.214
---
±0.083
| 0.243
---
±0.112
- |  | 0.227
---
±0.141
- |  | 0.217
---
±0.171
| 0.199
---
±0.080
| $\mathrm{RMSE}_{\mathbf{b}^{g}}$
---
[$rad/s$]
- |  | 0.031
---
±0.025
- |  | 0.048
---
±0.100
| 0.049
---
±0.012
- |  | 0.086
---
±0.168
- |  | 0.098
---
±0.209
| 0.103
---
±0.073
| $\theta_{\textbf{g}}$
---
[$rad$]
- |  | 0.372
---
±0.308
- |  | 0.299
---
±0.236
| 0.168
---
±0.089
- |  | 0.273
---
±0.204
- |  | 0.275
---
±0.215
| 0.167
---
±0.086

Table 1: Error Metrics on ICL- NUIM for different operating modes. (DIFODO∗)
stands for our implementation of the method in [7]. (∗∗) indicates the output
of our DIFODO∗ implementation between pairs of frames.

Figure 6: Motion estimation in an office scene from the ICL-NUIM dataset in
two different times $t_{1}$ and $t_{2}$. (a) 3D representation of the scene.
(b) Motion estimation of the objects in the scene. Every velocity is
represented by a red arrow on each point. (c) Zoomed-in areas.

| 2-frames | 3-frames | 4-frames | 5-frames
---|---|---|---|---
|  | RGB-D
---
(∗∗)
| RGB-D-I
---
(ours)
| RGB-D
---
(DIFODO∗)
| RGB-D-I
---
(ours)
| RGB-D
---
(∗∗)
| RGB-D-I
---
(ours)
| RGB-D
---
(∗∗)
| RGB-D-I
---
(ours)
| $\mathrm{RMSE}_{{\mathbf{v}}}$
---
[$cm/s$]
| 36.993
---
±11.295
| 34.008
---
±8.438
| 22.453
---
±9.029
| 22.054
---
±8.766
| 22.453
---
±9.029
| 21.808
---
±8.646
| 22.453
---
±9.029
| 21.924
---
±8.769
| $\mathrm{RMSE}_{{\boldsymbol{\omega}}}$
---
[$rad/s$]
| 0.165
---
±0.090
| 0.184
---
±0.102
| 0.172
---
±0.091
| 0.169
---
±0.092
| 0.172
---
±0.091
| 0.168
---
±0.092
| 0.172
---
±0.091
| 0.167
---
±0.092
| $\mathrm{RMSE}_{\mathbf{b}^{a}}$
---
[$cm/s^{2}$]
- |  | 0.121
---
±0.001
- |  | 0.102
---
±0.028
- |  | 0.097
---
±0.014
- |  | 0.113
---
±0.024
| $\mathrm{RMSE}_{\mathbf{b}^{g}}$
---
[$rad/s$]
- |  | 0.015
---
±0.001
- |  | 0.015
---
±0.001
- |  | 0.015
---
±0.001
- |  | 0.015
---
±0.001
| $\theta_{\textbf{g}}$
---
[$rad$]
- |  | 1.550
---
±0.440
- |  | 0.150
---
±0.104
- |  | 0.115
---
±0.072
- |  | 0.095
---
±0.061

Table 2: Error Metrics on OpenLORIS-Scene for different operating modes.
(DIFODO∗) stands for our implementation of the method in [7]. (∗∗) indicates
the output of our DIFODO∗ implementation between pairs of frames.

##  6 Experiments

###  6.1 Setup

Public benchmarks that provide IMU, color and depth images are scarce. We
chose to evaluate our proposal on an extended version of the living room
sequences in the ICL-NUIM dataset [5]. ICL-NUIM is a synthetic photorealistic
dataset that provides ground truth poses as well as 3D scene models to
benchmark reconstruction and/or localization approaches. As ICL-NUIM does not
provide IMU data, in a manner similar to [11], we fit splines to the ground
truth poses to simulate continuous trajectories and simulated IMU measurements
from them. We use the IMU model described in [18] and the same IMU parameters
as [12]. We also evaluated our RGB-D-inertial flow in the OpenLORIS-Scene
datasets [28], in which data are collected in real-world indoor scenes, for
multiple times in each place to include natural scene changes in everyday
scenarios. RGB-D images and IMU measurements from a RealSense D435i are
provided. The ground truth trajectory was recorded by an OptiTrack MCS, that
tracked artificial markers deployed on the Segway robot used to record the
data.

As metrics, we use the Root-Mean-Square-Error ($\mathrm{RMSE}$) for the
velocities $\mathbf{v}$ and ${\boldsymbol{\omega}}$ and biases
$\mathbf{b}^{a}$ and $\mathbf{b}^{g}$. For the gravity vector, we use the
angle
$\theta_{\textbf{g}}=\mathrm{cos}^{-1}\left(\nicefrac{{\hat{\textbf{g}}\cdot\textbf{g}_{gt}}}{{\|\hat{\textbf{g}}\|\|\textbf{g}_{gt}\|}}\right)$
between the ground truth and estimated gravity direction.

###  6.2 Results

ICL-NUIM. In this experiment we use the living room sequences in the ICL-NUIM
dataset, and we run our RGB-D-I flow based method against the so-called DIFODO
[7], based on RGB-D flow. Note that, as [7] does not provide code, we used our
own implementation based on the description in the paper. We run both scene
flow methods for different estimation modes: 2-frames, 3-frames, 4-frames and
5-frames. Marginalization is not done here. The results are obtained over the
entire dataset, taking as starting frame one out of every 2 which gives us
more than 400 subsequences. In this experiment we have made 10 runs on the
complete dataset. The results (specifically, the mean ± the standard deviation
of the RMSE for the estimated states) are shown in Table 1, in which the best
result per estimation mode is boldfaced. Note that using inertial measurements
improves the accuracy in both linear and angular velocity estimates. We can
also observe that the errors are higher for the 2-frames case than for the
rest, which shows how the information from additional frames is leveraged to
estimate the inertial states. Note also how the standard deviation of the
errors is reduced when inertial sensing is used, which indicates a higher
robustness in challenging cases.

ICL-NUIM + Marginalization. In this experiment, we show the effect of the
marginalization in the living room sequences of ICL-NUIM, using two different
sliding windows: 3 frames and 5 frames. In both cases we marginalize one frame
only. The results are shown in the columns of Table 1 titled as RGB-D-I
(ours+M). It can be seen how the gravity vector estimation is improved in
comparison to the case without marginalization. Also, a small improvement
occurs in the linear and angular velocities errors.

Fig. 6 shows qualitative results for motion estimation. The office scene
consists on a computer on a desk and luminaires. Fig. 6(a) shows a point cloud
extracted from the RGB-D data, and Fig. 6(b) displays the scene flow (red
arrows represent the velocity in each point). For better appreciation, Fig.
6(c) zooms in some areas. Observe that our approach estimates a smooth flow
even in textureless areas such as the background wall.

OpenLORIS-Scene. In this experiment we consider the office-1 sequence in the
OpenLORIS-Scene Dataset, where the robot moves along a U-shape route. As in
previous experiments, we compare our RGB-D-I flow-based motion estimation
against RGB-D-only motion estimation in 4 different optimization modes:
2-frames, 3-frames, 4-frames and 5-frames, all of them without
marginalization. The experiment is performed over the entire dataset, using
one frame out of every five as starting point, which gave us more than 400
experiments. Table 2 shows the mean ± the standard deviation of the RMSE for
the results of such experiments. As in the synthetic case, it can be observed
that using inertial measurements improves the estimation results in both
linear and angular velocities. It can also be observed how, as the
optimization window grows, the errors of the inertial states are also smaller.

##  7 Conclusions

In this work we present a novel camera motion estimation based on RGB-D-I
scene flow. Specifically, we formulate the fusion of RGB-D and inertial data
as a joint optimization using scene flow residuals and pre-integrated IMU
residuals, weighted by their corresponding covariances. We also consider the
marginalization of old states in order to keep a compact optimization. We
evaluated our approach on a synthetic dataset, ICL-NUIM, and on a real
dataset, OpenLORIS, both publicly available. Our results quantify the
improvement that inertial fusion can offer to RGB-D scene flow techniques.

## References

  * Armesto et al. [2007] Leopoldo Armesto, Josep Tornero, and Markus Vincze.  Fast ego-motion estimation with multi-rate fusion of inertial and vision.  _The International Journal of Robotics Research_ , 26(6):577–589, 2007.
  * Beauchemin and Barron [1995] S. S. Beauchemin and J. L. Barron.  The computation of optical flow.  _ACM Computing Surveys_ , 27(3):433–466, 1995.
  * Campos et al. [2021] Carlos Campos, Richard Elvira, Juan J Gómez Rodríguez, José MM Montiel, and Juan D Tardós.  Orb-slam3: An accurate open-source library for visual, visual–inertial, and multimap slam.  _IEEE Transactions on Robotics_ , 37(6):1874–1890, 2021.
  * Forster et al. [2017] Christian Forster, Luca Carlone, Frank Dellaert, and Davide Scaramuzza.  On-Manifold Preintegration for Real-Time Visual–Inertial Odometry.  _IEEE Transactions on Robotics_ , 33(1):1–21, 2017.
  * Handa et al. [2014] Ankur Handa, Thomas Whelan, John McDonald, and Andrew J. Davison.  A benchmark for rgb-d visual odometry, 3d reconstruction and slam.  In _2014 IEEE International Conference on Robotics and Automation (ICRA)_ , pages 1524–1531, 2014.
  * Hertzberg et al. [2013] Christoph Hertzberg, René Wagner, Udo Frese, and Lutz Schröder.  Integrating generic sensor fusion algorithms with sound state representations through encapsulation of manifolds.  _Information Fusion_ , 14(1):57–77, 2013.
  * Jaimez and Gonzalez-Jimenez [2015] Mariano Jaimez and Javier Gonzalez-Jimenez.  Fast visual odometry for 3-d range sensors.  _IEEE Transactions on Robotics_ , 31(4):809–822, 2015.
  * Jaimez et al. [2015] Mariano Jaimez, Mohamed Souiai, Javier Gonzalez-Jimenez, and Daniel Cremers.  A primal-dual framework for real-time dense rgb-d scene flow.  In _2015 IEEE international conference on robotics and automation (ICRA)_ , pages 98–104. IEEE, 2015.
  * Jaimez et al. [2018] Mariano Jaimez, Javier Monroy, Manuel Lopez-Antequera, and Javier Gonzalez-Jimenez.  Robust planar odometry based on symmetric range flow and multiscan alignment.  _IEEE Transactions on Robotics_ , 34(6):1623–1635, 2018.
  * Kerl et al. [2013] Christian Kerl, Jurgen Sturm, and Daniel Cremers.  Robust odometry estimation for RGB-D cameras.  In _2013 IEEE International Conference on Robotics and Automation_ , pages 3748–3754, Karlsruhe, Germany, 2013. IEEE.
  * Kerl et al. [2015] Christian Kerl, Jörg Stückler, and Daniel Cremers.  Dense continuous-time tracking and mapping with rolling shutter rgb-d cameras.  In _2015 IEEE International Conference on Computer Vision (ICCV)_ , pages 2264–2272, 2015.
  * Laidlow et al. [2017a] Tristan Laidlow, Michael Bloesch, Wenbin Li, and Stefan Leutenegger.  Dense rgb-d-inertial slam with map deformations.  In _2017 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS)_ , pages 6741–6748, 2017a.
  * Laidlow et al. [2017b] Tristan Laidlow, Michael Bloesch, Wenbin Li, and Stefan Leutenegger.  Dense RGB-D-inertial SLAM with map deformations.  In _2017 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS)_ , pages 6741–6748, Vancouver, BC, 2017b. IEEE.
  * Lee et al. [2020] Keuntaek Lee, Jason Gibson, and Evangelos A Theodorou.  Aggressive perception-aware navigation using deep optical flow dynamics and pixelmpc.  _IEEE Robotics and Automation Letters_ , 5(2):1207–1214, 2020.
  * Leutenegger et al. [2015] Stefan Leutenegger, Simon Lynen, Michael Bosse, Roland Siegwart, and Paul Furgale.  Keyframe-based visual–inertial odometry using nonlinear optimization.  _The International Journal of Robotics Research_ , 34(3):314–334, 2015.
  * Li et al. [2020] Xiuxiu Li, Yanjuan Liu, Haiyan Jin, Lei Cai, and Jiangbin Zheng.  RGBD Scene Flow Estimation with Global Nonrigid and Local Rigid Assumption.  _Discrete Dynamics in Nature and Society_ , 2020:1–9, 2020.
  * Nießner et al. [2014] Matthias Nießner, Angela Dai, and Matthew Fisher.  Combining inertial navigation and icp for real-time 3d surface reconstruction.  In _Eurographics (Short Papers)_ , pages 13–16. Citeseer, 2014.
  * Nikolic et al. [2016] Janosch Nikolic, Paul Furgale, Amir Melzer, and Roland Siegwart.  Maximum likelihood identification of inertial sensor noise model parameters.  _IEEE Sensors Journal_ , 16(1):163–176, 2016.
  * Nitsche et al. [2020] Matas Nitsche, Facundo Pessacg, and Javier Civera.  Visual-inertial teach and repeat.  _Robotics and Autonomous Systems_ , 131:103577, 2020.
  * Pire et al. [2017] Taihú Pire, Thomas Fischer, Gastón Castro, Pablo De Cristóforis, Javier Civera, and Julio Jacobo Berlles.  S-ptam: Stereo parallel tracking and mapping.  _Robotics and Autonomous Systems_ , 93:27–42, 2017.
  * Porzi et al. [2020] Lorenzo Porzi, Markus Hofinger, Idoia Ruiz, Joan Serrat, Samuel Rota Bulo, and Peter Kontschieder.  Learning multi-object tracking and segmentation from automatic annotations.  In _Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition_ , pages 6846–6855, 2020.
  * Prest et al. [2012] Alessandro Prest, Christian Leistner, Javier Civera, Cordelia Schmid, and Vittorio Ferrari.  Learning object class detectors from weakly annotated video.  In _2012 IEEE Conference on computer vision and pattern recognition_ , pages 3282–3289. IEEE, 2012.
  * Quiroga et al. [2014] Julian Quiroga, Thomas Brox, Frédéric Devernay, and James Crowley.  Dense Semi-rigid Scene Flow Estimation from RGBD Images.  In _Computer Vision – ECCV 2014_ , pages 567–582. Springer International Publishing, Cham, 2014.  Series Title: Lecture Notes in Computer Science.
  * Rosinol et al. [2023] Antoni Rosinol, John J Leonard, and Luca Carlone.  Nerf-slam: Real-time dense monocular slam with neural radiance fields.  In _2023 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS)_. IEEE, 2023.
  * Schnider et al. [2023] Yannick Schnider, Stanisław Woźniak, Mathias Gehrig, Jules Lecomte, Axel Von Arnim, Luca Benini, Davide Scaramuzza, and Angeliki Pantazi.  Neuromorphic optical flow and real-time implementation with event cameras.  In _Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition_ , pages 4128–4137, 2023.
  * Shan et al. [2020] Tixiao Shan, Brendan Englot, Drew Meyers, Wei Wang, Carlo Ratti, and Daniela Rus.  Lio-sam: Tightly-coupled lidar inertial odometry via smoothing and mapping.  In _2020 IEEE/RSJ international conference on intelligent robots and systems (IROS)_ , pages 5135–5142. IEEE, 2020.
  * Shan et al. [2019] Zeyong Shan, Ruijian Li, and Sören Schwertfeger.  RGBD-Inertial Trajectory Estimation and Mapping for Ground Robots.  _Sensors_ , 19(10):2251, 2019.
  * Shi et al. [2020] Xuesong Shi, Dongjiang Li, Pengpeng Zhao, Qinbin Tian, Yuxin Tian, Qiwei Long, Chunhao Zhu, Jingwei Song, Fei Qiao, Le Song, Yangquan Guo, Zhigang Wang, Yimin Zhang, Baoxing Qin, Wei Yang, Fangshi Wang, Rosa H. M. Chan, and Qi She.  Are we ready for service robots? the OpenLORIS-Scene datasets for lifelong SLAM.  In _2020 International Conference on Robotics and Automation (ICRA)_ , pages 3139–3145, 2020.
  * Spies et al. [2002] Hagen Spies, Bernd Jähne, and John L. Barron.  Range flow estimation.  _Comput. Vis. Image Underst._ , 85(3):209–231, 2002.
  * Teed and Deng [2021] Zachary Teed and Jia Deng.  Droid-slam: Deep visual slam for monocular, stereo, and rgb-d cameras.  _Advances in neural information processing systems_ , 34:16558–16569, 2021.
  * Vedula et al. [1999] S. Vedula, S. Baker, P. Rander, R. Collins, and T. Kanade.  Three-dimensional scene flow.  In _Proceedings of the Seventh IEEE International Conference on Computer Vision_ , pages 722–729 vol.2, Kerkyra, Greece, 1999. IEEE.
  * Zhai et al. [2021] Mingliang Zhai, Xuezhi Xiang, Ning Lv, and Xiangdong Kong.  Optical flow and scene flow estimation: A survey.  _Pattern Recognition_ , 114:107861, 2021.
  * Zhang et al. [2020] Tianwei Zhang, Huayan Zhang, Yang Li, Yoshihiko Nakamura, and Lei Zhang.  Flowfusion: Dynamic dense rgb-d slam based on optical flow.  In _2020 IEEE International Conference on Robotics and Automation (ICRA)_ , pages 7322–7328. IEEE, 2020.
  * Zihao Zhu et al. [2018] Alex Zihao Zhu, Liangzhe Yuan, Kenneth Chaney, and Kostas Daniilidis.  Unsupervised event-based optical flow using motion compensation.  In _Proceedings of the European Conference on Computer Vision (ECCV) Workshops_ , 2018.