Simultaneous localization and mapping for autonomous bicycles

Abstract

Bicycles and similar two-wheeled vehicles are a scientifically interesting category of means of transportation, whose simultaneous localization and mapping problem, in contrast to other ground vehicles, has yet to be examined in depth and fully resolved. In this article, we introduce for the first time a comprehensive theoretical framework for the application of simultaneous localization and mapping specifically for autonomous bicycles and akin vehicles, based on the kinematics and dynamics models of the bicycle and the specific effects they have on the motion/odometry and measurement models that are essential for the solution of the simultaneous localization and mapping problem. In addition, we present our laboratory’s first autonomous bicycle platform and its functional sensor system and sensor rotation pattern, specifically adjusted to the special characteristics of bicycles. Moreover, we investigate the effect of the bicycle frame roll, the main uncertainty factor of the bicycle simultaneous localization and mapping problem, on the overall simultaneous localization and mapping performance. The experimental results performed on our bicycle platform verify the potency of our proposed modeling and simultaneous localization and mapping application framework and provide further insight on future improvements for the two-wheeled vehicle simultaneous localization and mapping problem.

Introduction

Simultaneous localization and mapping (SLAM) is a fairly well-known problem for autonomous robots that still presents numerous difficulties for researchers, even though a variety of efficient solutions have been proposed mainly after the dawn of the 21st century. SLAM applications are widespread, and its benefits have been clearly demonstrated for the mapping and navigation of autonomous robots. Although a great deal of effort has focused in the recent years in fully resolving the problem of driverless cars, other categories of autonomous vehicles and their specific requirements for a smooth and efficient application of SLAM have not been sufficiently examined. Such a type of vehicles is the autonomous bicycle or motorcycle. Compared to unmanned cars, bicycles are more flexible and can freely pass through narrow spaces such as alleys and country roads, while they constitute a widely used means of transportation, mainly due to their low cost and ease of use. In contrast to four-wheeled vehicles, the virtues of autonomous bicycles and similar two-wheeled vehicles in general have yet to be thoroughly presented and understood. Lightweight, easily controlled, fast in changing direction, cost-, space-, and energy-efficient bicycles have a plethora of advantages and possible applications once they achieve autonomous behavior. In order to reach that stage, however, the achievement of reliable SLAM, that is, constructing a model of their surroundings and navigating with stability, speed, and safety within the map based on that model, is a necessity.

That necessity can also be understood, if we consider the following. Firstly, SLAM can handle difficult and highly dynamic environment conditions and sensor fusion techniques under more challenging setups, for example, no GPS, low-quality sensors, than examined by previous localization and mapping methods, which is usual in the case of more simple vehicle configurations that befit bicycles. Secondly, by finding loop closures with SLAM, autonomous bicycles can understand the real environment topology and are able to find shortcuts between different map locations, while the metric map reconstruction can prevent erroneous loop closures from happening, enhancing robustness. Finally, due to numerous autonomous bicycle applications being in need of a globally consistent map, such as goods deliveries or intra-city navigation, the implementation of bicycle SLAM is deemed crucial. More specifically, a successful SLAM should account for the two-wheeled vehicles’ complex kinematics and dynamics, include a sensor system that can use the bicycles’ one wheel odometry, account for the constant vibrations of the bicycle frame, especially on uneven types of terrain, and make sure that real-time performance can be achieved for navigation and obstacle avoidance purposes.

Regardless of the SLAM methodology applied, there are two vital components for their application based on the Bayesian recursion model. ¹ First, the reliable prediction or registration of the robot’s motion, given an accurate motion or odometry model, respectively, in order to execute the “sample” phase. And second, the collection of measurements of the surrounding environment, using a variety of sensors with known measurement models, in order to execute the “update” phase. The theoretical basis of autonomous vehicles and obstacle avoidance has been firmly set over the past decades, mainly following the principles of automatic control. ² Publications of the last decade have introduced several solutions for the autonomous cars’ SLAM problem. ^3,4 Various teams proposed techniques for environment perception of autonomous vehicles using range-bearing finding sensors, ^5,6 while others analyzed the solution to the visual SLAM problem. ^7,8 Sensor fusion techniques have also been implemented by numerous researchers for outdoor environments, ⁹ whereas RGB-D sensors have not yet been fully utilized due to their still poor results in outdoor environments.

One of the first efforts pertaining specifically to autonomous bicycles was ¹⁰ examining the stabilization and control methodologies of such two-wheeled vehicles. Other teams ^11,12 have studied bicycles as well, but mainly from the control aspect. Yi et al. provided a thorough analysis of the bicycle’s kinematics and dynamics models, ^13,14 while an applicable version of these models for a real-world implementation of a self-stabilizing bicycle was presented in the work by He et al., ¹⁵ and obstacle avoidance methodologies were examined by Stasinopoulos et al. ¹⁶ Further real-world applications of autonomous bicycle perception are relatively scarce, while no team has analyzed the two-wheeled SLAM problem before, mainly due to the lack of a sound theoretical framework describing their perception and SLAM aspects. A simplistic approach could suggest that the same principles that apply for four-wheeled vehicles, as seen in research efforts by teams for both two-dimensional ^17,18 or three-dimensional ^19,20 environment representations, could be used to analyze the SLAM problem for autonomous bicycles. However, in reality, a significant number of fundamental differences between these two types of vehicles exist that demand the redefinition and adaptation of the problem to conform to the special characteristics of these autonomous vehicles.

This article provides a thorough theoretical analysis of the unique characteristics of autonomous bicycles that differentiate them from other mobile robots with regards to the application of SLAM. Our main contributions comprise the modification of the differential motion model for mobile robots based on the complex bicycle kinematics and dynamics, two types of odometry models along with the uncertainty factors they introduce, the introduction of an applicable sensor design and rotation pattern for onboard perception, and the specially modified measurement model for landmark detection. Our SLAM application solution can be extended for use by all similar two-wheeled vehicles. In addition, we present for the first time our laboratory’s autonomous platform prototype and its fully functional control, odometry registration and perception systems, capable of SLAM application. Furthermore, we analyze how the rate of change of the bicycle frame roll angle and the overall SLAM performance are connected and verify this hypothesis through our experiments.

The remainder of the article is organized as follows: “Bicycle model” section introduces the kinematics and dynamics models for the bicycle, while “Motion and odometry models” and “Environment sensing system and measurement model” sections present our proposed motion/odometry model and sensory system/measurement model, respectively. “Implications on SLAM application” section provides discussion on the main difficulties of the two-wheeled vehicle SLAM problem and establishes a connection between the bicycle roll change rate and the SLAM algorithm performance. “Autonomous bicycle prototype and experimental setting” and “Experimental results” sections describe our autonomous bicycle prototype and the experimental system and results that verify our bicycle SLAM problem solution’s potency and our hypothesis regarding performance. Finally, “Conclusion” section presents the conclusions of our effort and possible further extensions.

Bicycle model

As mentioned above, one of the most accurate models was presented in the work of Yi et al., ^13,14 where the bicycle is considered as a two-part platform, a rear frame, and a steering mechanism. Figure 1 shows a schematic of the vehicle. They assume that (1) the wheel ground is a point contact and thickness and geometry of the tire are neglected; (2) the bicycle body frame is considered a point mass; and (3) the bicycle moves on a flat plane with no vertical motion.

Figure 1.

A kinematic diagram of the bicycle dynamics. (a) Schematic view of a bicycle on a plane and (b) top-view of the bicycle, as presented in the work of Yi et al. ¹⁴

In this model, X, Y, Z is the ground-fixed coordinate system, x, y, z the wheelbase line moving coordinate system, $x_w,y_w,z_w$ the front wheel plane coordinate system, and $x_B,y_B,z_B$ the rear frame body coordinate system. Let C₁, C₂ be the front and rear wheel/ground contact points, and F_ix, F_iy, F_iz the (front and rear) wheel contact forces in the x, y, z-axis directions, where i = f, r. Also, let υ_f, υ_r be the velocity vectors of the front and rear wheel contact points, respectively, and υ_ix, υ_iy the (front and rear) wheel contact point velocities along the x- and y-axis directions, respectively, with i = f, r. Moreover, υ_fxw, υ_fyw are the front wheel contact point C₁ velocities along the front wheel directions, υ_X, υ_Y are the rear wheel contact point C₂ velocities along the navigation directions, respectively, and ω_f, ω_r are the wheel angular velocities of the front and rear wheels, respectively. In addition, υ_G is the velocity vector of the bicycle frame (with rear wheel set), while λ_i, γ_i are the (front and rear) tire slip ratio and angle, respectively, where i = f, r. Then, ϕ, ψ are the rear frame roll and yaw angles, respectively, and φ, φ_g are the steering and projected steering angles, respectively. Furthermore, σ is the kinematic steering angle variable, m the mass of the rear frame and wheel, l the bicycle wheel base, h the height of the bicycle center of mass, and r the (front and rear) wheel radius. Finally, if ξ is the steering axis caster angle and R the radius of the trajectory of point C₂, respectively, then the bicycle motion and dynamics can be described as follows.

The kinematic steering variable σ is defined as $\sigma =\text{tan} {\phi }_g =\frac{l}{R}$ . From the geometry of the front wheel steering mechanism, the following relationship is found $\text{tan}{\phi }_g{c}_{\varphi }=\text{tan}\phi c_{\xi }$ , where $c_x:=\cos x$ , $s_x:=\text{sin}x$ for angle x. If small roll and steering angles are assumed, then in approximation $\dot{\sigma }{c}_{\varphi }=\dot{\phi }{c}_{\xi }$ . From the top-view, Figure 1(b), the relationship between (x, y) and yaw angle ψ is calculated as

$\left[\begin{array}{c}\dot{x}\\ \dot{y}\end{array}\right]=\left[\begin{array}{cc}\cos \psi & -\text{sin}\psi \\ \text{sin}\psi & \cos \psi \end{array}\right]\left[\begin{array}{c}{\upsilon }_{Gx}\\ {\upsilon }_{Gy}\end{array}\right]$ 1

The nonholonomic constraint of the rear wheel implies that

${\upsilon }_{rx}=R\dot{\psi }=\frac{l}{\sigma }\dot{\psi }$ 2

In such a way, the mass center velocity is derived as

${\upsilon }_G=\left({\upsilon }_{rx}-h\dot{\psi }{s}_{\varphi }\right)i+\left({\upsilon }_{ry}+l\dot{\psi }+r\dot{\varphi }{s}_{\varphi }\right)j+h\dot{\varphi }{s}_{\varphi }k$ 3

In a similar fashion, by introducing the Lagrangian L and the constrained Lagrangian L_c, the bicycle’s roll, longitudinal, and lateral dynamics equations are calculated. Let $\dot{q}:= {\left[\dot{\varphi } {\upsilon }_{rx} {\upsilon }_{ry}\right]}^{\text{T}}$ denote the generalized velocity of the bicycle and thus, the dynamics equations can be written in a compact matrix form as

$\begin{array}{c}\dot{\sigma } = {\omega }_{\sigma }\\ M\left[q\right] = K_m + B_m\left[\begin{array}{c}{\omega }_{\sigma }\\ F_{fx}\\ F_{fy}\\ F_{rx}\\ F_{ry}\end{array}\right]\end{array}$ 4

where ω_σ is the virtual steering angular velocity input, and the matrices

$\begin{array}{l}M= \left[\begin{array}{ccc}{h}^2& \frac{bh\sigma }{l}{c}_{\varphi }& h{c}_{\varphi }\\ \frac{bh\sigma }{l}{c}_{\varphi }& {\left(1-\frac{h\sigma }{l}{s}_{\varphi }\right)}^2+\frac{b^2{\sigma }^2}{l^2}& \frac{b\sigma }{l}\\ h{c}_{\varphi }& \frac{b\sigma }{l}& 1\end{array}\right],K_m= \left[\begin{array}{c}-\left(1-\frac{h\sigma }{l}{s}_{\varphi }\right)\frac{h\sigma c_{\varphi }}{l}{\upsilon }_{rx}{ }^2+g\left(h{s}_{\varphi }+\frac{l_{t}b{c}_{\xi }\sigma }{l}{c}_{\varphi }\right)\\ 2\left(1-\frac{h\sigma }{l}{s}_{\varphi }\right)\frac{h\sigma c_{\varphi }}{l}\dot{\varphi }{\upsilon }_{rx}+\frac{bh\sigma }{l}{s}_{\varphi }{\dot{\varphi }}^2-\frac{1}{m}{C}_d{\upsilon }_{rx}{ }^2\\ h{s}_{\varphi }{\dot{\varphi }}^2\end{array}\right],\\ \text{and}\hspace{0.17em}{B}_m=\left[\begin{array}{ccccc}-\frac{bh}{l}{c}_{\varphi }{\upsilon }_{rx}& 0& 0& 0& 0\\ A_m& -\frac{1}{m\sqrt{1+{\sigma }^2}}& -\frac{\sigma }{m\sqrt{1+{\sigma }^2}}& \frac{1}{m}& 0\\ -\frac{b}{l}{\upsilon }_{rx}& -\frac{\sigma }{m\sqrt{1+{\sigma }^2}}& -\frac{1}{m\sqrt{1+{\sigma }^2}}& 0& -\frac{1}{m}\end{array}\right]\hspace{0.17em}\text{with}\hspace{0.17em}{A}_m=2\left[\left(1-\frac{h\sigma }{l}{s}_{\varphi }\right)\frac{h}{l}{s}_{\varphi }-\frac{b^2{\sigma }^{ }}{l^2}\right]{\upsilon }_{rx}-\frac{b}{l}{\upsilon }_{ry}-\frac{bh{c}_{\varphi }}{l}\dot{\varphi }\end{array}$

It is clear that the control inputs in equation (4) are the virtual steering velocity ω_σ and the wheel traction/braking forces $F_{fx},F_{rx},F_{fy}$ , and F_ry. In the bicycle model presented above, the front and rear tires’ longitudinal slip ratio λ_s and the lateral slip ratio λ_γ are disregarded, which as explained in the study by Yi et al. ¹⁴ can be expressed as

${\lambda }_s= \frac{{\upsilon }_{cx}-r {\omega }_w}{{\upsilon }_{cx}} , {\lambda }_{\gamma }=\text{tan}\gamma =- \frac{{\upsilon }_{cy}}{{\upsilon }_{cx}}$ 5

where ω_w is the wheel angular velocity and γ is the sideslip angle, while υ_cx and υ_cy are the velocities on the x and y axes, respectively, of the contact points C₁ and C₂ for the front and rear wheels, respectively. In various occasions, however, this slip cannot be neglected for both the front and the rear wheel and thus, the front wheel longitudinal slip ratio and the rear wheel sideslip ratio are, respectively

${\lambda }_{fs}=1-\frac{r{\upsilon }_f}{{\upsilon }_{fxw}}=1-\frac{r\sqrt{1+{\sigma }^2}}{\left(1+{\sigma }^2\right){\upsilon }_{rx}+\sigma {\upsilon }_{ry}}{\omega }_f, {\lambda }_{r\gamma }=-\frac{{\upsilon }_{fy}}{{\upsilon }_{fx}}-\frac{r\text{tan}\xi c_{\varphi }{ }^2}{{\upsilon }_{rx}}{\omega }_{\sigma }+\sigma$ 6

Motion and odometry models

From the above description of the bicycle model and dynamics, it is obvious that the motion of a bicycle is fundamentally different from a four-wheeled vehicle’s. Considering the effect this model has on the registration of the bicycle’s displacement and change in orientation at every instant, that is, the registration of its odometry, we can observe from equation (3) that in order to calculate the Cartesian displacement ${\left[d{x}_G,d{y}_G\right]}^{\text{T}}$ we need to have knowledge of the center of mass velocity υ_G and the yaw angle of the bicycle frame ψ.

Motion model

As can be seen from the bicycle model above, the motion model from pose $x_{t-1}={\left[x_G,y_G,\psi \right]}^{\text{T}}$ to pose $x_t={\left[{x\prime }_G,{y\prime }_G,\psi \prime \right]}^{\text{T}}$ after time Δt cannot be described by the simple representation that uses a combination of the longitudinal velocity υ_rx and the rotational velocity $\dot{\psi }$ , as proposed by Thrun et al. ¹ for generic vehicles

$\left[\begin{array}{c}{x\prime }_G\\ {y\prime }_G\\ \psi \prime \end{array}\right]=\left[\begin{array}{c}{x}_G\\ y_G\\ \psi \end{array}\right]+\left[\begin{array}{c}-\frac{ {\upsilon }_{rx}}{ \dot{\psi }}\text{sin}\psi +\frac{ {\upsilon }_{rx}}{ {\omega }_G}\text{sin}\left(\psi +\dot{\psi }\Delta t\right)\\ \frac{ {\upsilon }_{rx}}{ \dot{\psi }}\cos \psi -\frac{ {\upsilon }_{rx}}{ {\omega }_G}\cos \left(\psi +\dot{\psi }\Delta t\right)\\ \dot{\psi }\Delta t+\gamma \Delta t\end{array}\right]$ 7

where γ is the extra final rotation added to the model to prevent the degenerate behavior of the posterior distribution. That is because the motion is more complex in the case of the bicycle and it includes the lateral movement caused mainly by the roll angle. Therefore, the proper equation based on our model should be

$\left[\begin{array}{c}{x\prime }_G\\ {y\prime }_G\\ \psi \prime \end{array}\right]=\left[\begin{array}{c}{x}_G\\ y_G\\ \psi \end{array}\right]+\left[\begin{array}{c}\left({\upsilon }_{rx}-h\dot{\psi }{s}_{\varphi }\right)\Delta t\\ \left({\upsilon }_{ry}+l\dot{\psi }+r\dot{\varphi }{s}_{\varphi }\right)\Delta t\\ \dot{\psi }\Delta t\end{array}\right]$ 8

As can be observed, the main difference is that our motion also depends greatly on the roll angular velocity $\dot{\varphi }$ . Since these velocities can be calculated through the dynamics model presented by equation (4), in order to calculate the probability $p\text{(}{x}_t|u_t,x_{t-1})$ , where $u_t={\left[{\omega }_{\sigma },F_{fx},F_{fy},F_{rx},F_{ry}\right]}_t^T$ is the control input at time t, we need to define the errors between the commanded and the actual values of the control variables. That is because the robot motion is subject to noise and the actual virtual steering and traction/breaking forces on the wheels differ from the commanded ones, as seen below

$\left[\begin{array}{c}{\omega }_{\sigma \mathrm{,}\text{actual}}\\ F_{fx\mathrm{,}\text{actual}}\\ F_{fy\mathrm{,}\text{actual}}\\ F_{rx\mathrm{,}\text{actual}}\\ F_{ry\mathrm{,}\text{actual}}\end{array}\right]=\left[\begin{array}{c}{\omega }_{\sigma }\\ F_{fx}\\ F_{fy}\\ F_{rx}\\ F_{ry}\end{array}\right]+\left[\begin{array}{c}{\epsilon }_{{\alpha }_1\left|{\omega }_{\sigma }\right|+{\alpha }_2\left|F_{fx}\right|+{\alpha }_3\left|F_{fy}\right|+{\alpha }_4\left|F_{rx}\right|+{\alpha }_5\left|F_{ry}\right|}\\ {\epsilon }_{{\alpha }_6\left|{\omega }_{\sigma }\right|+{\alpha }_7\left|F_{fx}\right|+{\alpha }_8\left|F_{fy}\right|+{\alpha }_9\left|F_{rx}\right|+{\alpha }_{10}\left|F_{ry}\right|}\\ {\epsilon }_{{\alpha }_{11}\left|{\omega }_{\sigma }\right|+{\alpha }_{12}\left|F_{fx}\right|+{\alpha }_{13}\left|F_{fy}\right|+{\alpha }_{14}\left|F_{rx}\right|+{\alpha }_{15}\left|F_{ry}\right|}\\ {\epsilon }_{{\alpha }_{16}\left|{\omega }_{\sigma }\right|+{\alpha }_{17}\left|F_{fx}\right|+{\alpha }_{18}\left|F_{fy}\right|+{\alpha }_{19}\left|F_{rx}\right|+{\alpha }_{20}\left|F_{ry}\right|}\\ {\epsilon }_{{\alpha }_{21}\left|{\omega }_{\sigma }\right|+{\alpha }_{22}\left|F_{fx}\right|+{\alpha }_{23}\left|F_{fy}\right|+{\alpha }_{24}\left|F_{rx}\right|+{\alpha }_{25}\left|F_{ry}\right|}\end{array}\right]$ 9

where ε_b is a zero-mean error variable with variance b, such as a normal distribution with zero mean μ = 0 and variance σ² = b. In our model, the variance of the error is proportional to the commanded control variable. The parameters α₁ to α₂₅ (with ${\alpha }_i\ge 0$ for $i=1,\dots \mathrm{,25}$ ) are bicycle-specific error parameters and they model the accuracy of the bicycle. The less accurate the steering and traction/breaking systems on a bicycle, the larger these parameters.

Under our error model, specified in equation (9), these motion errors have the probabilities shown in the set below

$\left\{\begin{array}{c}{\epsilon }_{{\alpha }_1\left|{\omega }_{\sigma }\right|+{\alpha }_2\left|F_{fx}\right|+{\alpha }_3\left|F_{fy}\right|+{\alpha }_4\left|F_{rx}\right|+{\alpha }_5\left|F_{ry}\right|}\cdot {\omega }_{\sigma \mathrm{,}\text{error}}\mathrm{,}\\ {\epsilon }_{{\alpha }_6\left|{\omega }_{\sigma }\right|+{\alpha }_7\left|F_{fx}\right|+{\alpha }_8\left|F_{fy}\right|+{\alpha }_9\left|F_{rx}\right|+{\alpha }_{10}\left|F_{ry}\right|}\cdot F_{fx\mathrm{,}\text{error}}\mathrm{,}\\ {\epsilon }_{{\alpha }_{11}\left|{\omega }_{\sigma }\right|+{\alpha }_{12}\left|F_{fx}\right|+{\alpha }_{13}\left|F_{fy}\right|+{\alpha }_{14}\left|F_{rx}\right|+{\alpha }_{15}\left|F_{ry}\right|}\cdot F_{fy\mathrm{,}\text{error}}\mathrm{,}\\ {\epsilon }_{{\alpha }_{16}\left|{\omega }_{\sigma }\right|+{\alpha }_{17}\left|F_{fx}\right|+{\alpha }_{18}\left|F_{fy}\right|+{\alpha }_{19}\left|F_{rx}\right|+{\alpha }_{20}\left|F_{ry}\right|}\cdot F_{rx\mathrm{,}\text{error}}\mathrm{,}\\ {\epsilon }_{{\alpha }_{21}\left|{\omega }_{\sigma }\right|+{\alpha }_{22}\left|F_{fx}\right|+{\alpha }_{23}\left|F_{fy}\right|+{\alpha }_{24}\left|F_{rx}\right|+{\alpha }_{25}\left|F_{ry}\right|}\cdot F_{ry\mathrm{,}\text{error}}\end{array}\right\}$ 10

Since we assume independence between the different sources of error, the desired probability $p\text{(}{x}_t|u_t,x_{t-1})$ is the product of these individual errors

$\begin{array}{ll}p\text{(}{x}_t|u_t,x_{t-1})=\hfill & \hspace{0.17em}\hspace{0.17em}{\epsilon }_{{\alpha }_1\left|{\omega }_{\sigma }\right|+{\alpha }_2\left|F_{fx}\right|+{\alpha }_3\left|F_{fy}\right|+{\alpha }_4\left|F_{rx}\right|+{\alpha }_5\left|F_{ry}\right|}\cdot {\omega }_{\sigma \mathrm{,}\text{error}}\cdot {\epsilon }_{{\alpha }_6\left|{\omega }_{\sigma }\right|+{\alpha }_7\left|F_{fx}\right|+{\alpha }_8\left|F_{fy}\right|+{\alpha }_9\left|F_{rx}\right|+{\alpha }_{10}\left|F_{ry}\right|}\cdot F_{fx\mathrm{,}\text{error}}\hfill \\ \hfill & \cdot \hspace{0.17em}{\epsilon }_{{\alpha }_{11}\left|{\omega }_{\sigma }\right|+{\alpha }_{12}\left|F_{fx}\right|+{\alpha }_{13}\left|F_{fy}\right|+{\alpha }_{14}\left|F_{rx}\right|+{\alpha }_{15}\left|F_{ry}\right|}\cdot F_{fy\mathrm{,}\text{error}}\cdot {\epsilon }_{{\alpha }_{16}\left|{\omega }_{\sigma }\right|+{\alpha }_{17}\left|F_{fx}\right|+{\alpha }_{18}\left|F_{fy}\right|+{\alpha }_{19}\left|F_{rx}\right|+{\alpha }_{20}\left|F_{ry}\right|}\cdot F_{rx\mathrm{,}\text{error}}\hfill \\ \hfill & \cdot \hspace{0.17em}{\epsilon }_{{\alpha }_{21}\left|{\omega }_{\sigma }\right|+{\alpha }_{22}\left|F_{fx}\right|+{\alpha }_{23}\left|F_{fy}\right|+{\alpha }_{24}\left|F_{rx}\right|+{\alpha }_{25}\left|F_{ry}\right|}\cdot F_{ry\mathrm{,}\text{error}}\hfill \end{array}$ 11

Odometry model

Regarding the model where only odometry is registered instead of using the motion model to predict the behavior of the bicycle, we can obtain the main two components of our motion, namely the rear wheel longitudinal velocity υ_rx and the bicycle’s orientation-yaw angle ψ, using the main two following methods.

Velocity-orientation registration

Firstly, we can measure the rear wheel longitudinal velocity υ_rx by placing an encoder on the rear wheel resulting in a straightforward and efficient way. However, as observed in numerous occasions in real-world bicycle applications, the rear wheel longitudinal slip ratio λ_s and the lateral slip ratio λ_γ, as presented in equation (5), are not negligible. This can depend on the tire condition, since worn-out tires tend to slip more often, and on the road surface type and condition, whether we have a very slippery surface, road imperfections, or other items on the road, such as bumps or oil spills, that may increase the slipping probability. Although heavier and more stable vehicles may be able to perform better under these conditions, the lightweight structure of bicycles, coupled with the thinness of their tires, makes it vital to account for the longitudinal and lateral slip and integrate a specific uncertainty factor related to them when measuring υ_rx and ψ. Therefore, the uncertainty of the longitudinal velocity measurement, represented by the error variance ${\sigma }^2{ }_{ {\upsilon }_{rx}}$ , that usually only depends on the accuracy of the encoder measuring the speed of the rear wheel, represented by the error variance ${\sigma }^2{ }_{\text{vel,enc }}$ , and the lateral velocity uncertainty ${\sigma }^2{ }_{{\upsilon }_{ry}}$ will be, respectively

${\sigma }^2{ }_{ {\upsilon }_{rx}}={\sigma }^2{ }_{\text{ vel,enc}}\cdot c_1\frac{1}{\left|{\lambda }_s\right|} , {\sigma }^2{ }_{ {\upsilon }_{ry}}=c_2\frac{1}{\left|{\lambda }_{\gamma }\right|}$ 12

where c₁, c₂ are bicycle-specific parameters. Thus, the odometry error covariance matrix for the bicycle, regardless of which method is used to produce it, has to be adjusted accordingly to include the above.

Additionally, most four-wheeled, three-wheeled, or even two-wheeled—with the two wheels on the same axis—vehicles obtain their odometry mainly by means of attaching encoders to the differential drive system, that is, comparing the different rotation of the two wheels at either side of a main axis and calculating the instantaneous orientation of the robot. However, this technique cannot be implemented for bicycles, for which the rate of change of the yaw angle of the bicycle frame ψ according to equation (2) is

$\dot{\psi }= \frac{\sigma }{l} {\upsilon }_{rx}=\frac{\text{tan}{\phi }_g}{l}{\upsilon }_{rx}=\frac{\text{tan}\phi c_{\xi }}{l c_{\varphi }}{\upsilon }_{rx}\iff \Delta \psi =\frac{\text{tan}\phi c_{\xi }}{l c_{\varphi }}{\upsilon }_{rx}\Delta t$ 13

Therefore, since the caster angle η and the bicycle wheelbase length L are constants for each bicycle, the change in orientation of the bicycle during its movement and thus, the calculation of the virtual steering ω_σ, where ${\omega }_{\sigma }{c}_{\varphi }=\dot{\phi }{c}_{\xi }$ , can be approximated at each instant with (1) the usage of an encoder that registers the steering angle φ, (2) measurement of the longitudinal velocity of the rear wheel υ_r, as described above, and (3) measurement of the roll angle ϕ. The last measurement poses the greatest challenge, since it relates to the inherent characteristic of bicycle-like two-wheeled vehicles to roll while turning and is not a control variable. An accurate gyroscopic sensor can be used to register the roll angle ϕ at any instant, thus making the calculation of Δ_ψ possible. In that way, additional variances for the steering motor encoder ${\sigma }^2{ }_{\phi ,\text{enc} }$ and for the gyroscopic roll registration ${\sigma }^2{ }_{\varphi ,\text{gyro }}$ must be considered in order to calculate the final variance of the rotation measurement. If we approximate the motion of our bicycle based on the measurements by the rear wheel motor encoder, the steering motor encoder, and the gyroscopic sensor roll measurement, then from equation (1) we get

$\left[\begin{array}{c}{x\prime }_G\\ {y\prime }_G\\ \psi \prime \end{array}\right]=\left[\begin{array}{c}{x}_G\\ y_G\\ \psi \end{array}\right]+\left[\begin{array}{c}\left[\cos \left(\nu \right){\upsilon }_{rx}-\text{sin}\left(\nu \right){\upsilon }_{ry}\right]\Delta t\\ \left[\text{sin}\left(\nu \right){\upsilon }_{rx}+\cos \left(\nu \right){\upsilon }_{ry}\right]\Delta t\\ \frac{\text{tan}\phi c_{\xi }}{l c_{\varphi }}{\upsilon }_{rx}\Delta t\end{array}\right]$ 14

where $\nu =\left(\psi +\Delta \psi \right)=\psi +\frac{\text{tan}\phi c_{\xi }}{l c_{\varphi }}{\upsilon }_{rx}\Delta t$ , and we consider that the longitudinal rear wheel velocity υ_rx, the steering angle φ, and the roll angle ϕ remain constant for Δt. If we assume that the lateral rear wheel velocity can be represented by a zero-mean distribution of variance ${\sigma }^2{ }_{ {\upsilon }_{ry}}$ , then the measured variables can be defined as shown in the following set

$\left\{\begin{array}{c}{\upsilon }_{rx}={\upsilon }_{rx,\text{actual}}-{\epsilon }_{b_1\left|{\upsilon }_{rx}\right|+b_2\left|\phi \right|+b_3{\left|\varphi \right|}^{\prime }}\\ \phi ={\phi }_{\text{actual}}-{\epsilon }_{b_4\left|{\upsilon }_{rx}\right|+b_5\left|\phi \right|+b_6{\left|\varphi \right|}^{\prime }}\\ \varphi ={\varphi }_{\text{actual}}-{\epsilon }_{b_7\left|{\upsilon }_{rx}\right|+b_8\left|\phi \right|+b_9{\left|\varphi \right|}^{\prime }}\\ {\upsilon }_{ry}={\epsilon }_{{\sigma }^2{ }_{ {\upsilon }_{ry}}}\end{array}\right\}$ 15

while the errors due to equation (12) become

$\left\{\begin{array}{c}{\upsilon }_{rx}={\upsilon }_{rx,\text{actual}}-{\epsilon }_{d_1{\sigma }^2{ }_{ {\upsilon }_{rx}}\left|{\upsilon }_{rx}\right|+d_2{\sigma }^2{ }_{\phi ,\text{enc} }\left|\phi \right|+d_3{\sigma }^2{ }_{\varphi ,\text{gyro} }{\left|\varphi \right|}^{\prime }}\\ \phi ={\phi }_{\text{actual}}-{\epsilon }_{d_4{\sigma }^2{ }_{ {\upsilon }_{rx}}\left|{\upsilon }_{rx}\right|+d_5{\sigma }^2{ }_{\phi ,\text{enc} }\left|\phi \right|+d_6{\sigma }^2{ }_{\varphi ,\text{gyro} }{\left|\varphi \right|}^{\prime }}\\ \varphi ={\varphi }_{\text{actual}}-{\epsilon }_{d_7{\sigma }^2{ }_{ {\upsilon }_{rx}}\left|{\upsilon }_{rx}\right|+d_8{\sigma }^2{ }_{\phi ,\text{enc} }\left|\phi \right|+d_9{\sigma }^2{ }_{\varphi ,\text{gyro }}{\left|\varphi \right|}^{\prime }}\\ {\upsilon }_{ry}={\epsilon }_{d_{10}{\sigma }^2{ }_{ {\upsilon }_{ry}}}\end{array}\right\}$ 16

where ε_b are the error distributions as described above, and d₁ to d₁₀ are bicycle-specific error parameters. Therefore, if we assume the errors to be independent, the odometry model joint error probability is the product of the separate error probability distributions

$\begin{array}{ll}p\text{(}{x}_t|u_t,x_{t-1})=\hfill & {\epsilon }_{d_1{\sigma }^2{ }_{ {\upsilon }_{rx}}\left|{\upsilon }_{rx}\right|+d_2{\sigma }^2{ }_{\phi ,\text{enc} }|\phi |+d_3{\sigma }^2{ }_{\varphi ,\text{gyro} }|\varphi |}\cdot \left({\upsilon }_{rx,\text{actual}}-{\upsilon }_{rx}\right)\cdot {\epsilon }_{d_4{\sigma }^2{ }_{ {\upsilon }_{rx}}\left|{\upsilon }_{rx}\right|+d_5{\sigma }^2{ }_{\phi ,\text{enc} }|\phi |+d_6{\sigma }^2{ }_{\varphi ,\text{gyro} }|\varphi |}\cdot \left({\phi }_{\text{actual}}-\phi \right)\hfill \\ \hfill & \cdot {\epsilon }_{d_7{\sigma }^2{ }_{ {\upsilon }_{rx}}\left|{\upsilon }_{rx}\right|+d_8{\sigma }^2{ }_{\phi ,\text{enc} }|\phi |+d_9{\sigma }^2{ }_{\varphi ,\text{gyro} }|\varphi |}\cdot \left({\varphi }_{\text{actual}}-\varphi \right)\cdot {\epsilon }_{d_{10}{\sigma }^2{ }_{ {\upsilon }_{ry}}}\hfill \end{array}$ 17

Gyroscopic odometry

Of course, the above data can also be obtained using a gyroscopic sensor at the rear part of the bicycle that can measure its speed and orientation at each instant, thus providing a method of measuring the odometry that is uncorrelated with the bicycle and tire dynamics. As observed in real-world applications, however, the accuracy provided for the velocity and the yaw angle of the bicycle, that is, its orientation, are not as reliable as the measurements provided by the motor encoders described above. Therefore, measuring the velocity and orientation by means of a gyroscope can be used as an extra input to the system of an autonomous bicycle to validate the existing measurements. Otherwise, if we solely depend on the gyroscope, factors of gyroscope covariance for the velocity ${\upsilon }_{Gx,\text{gyro}}$ , ${\sigma }^2{ }_{\text{vel,gyro }}$ , and for the orientation ${\psi }_{\text{gyro}}$ , ${\sigma }^2{ }_{\psi ,\text{gyro} }$ , need to be integrated to our odometry error model as follows

$\left[\begin{array}{c}{x\prime }_G\\ {y\prime }_G\\ \psi \prime \end{array}\right]=\left[\begin{array}{c}{x}_G\\ y_G\\ \psi \end{array}\right]+\left[\begin{array}{c}\left[\cos \left(\psi +{\psi }_{\text{gyro}}\right){\upsilon }_{Gx,\text{gyro}}\right]t\\ \left[\text{sin}\left(\psi +{\psi }_{\text{gyro}}\right){\upsilon }_{Gx,\text{gyro}}\right]t\\ {\psi }_{\text{gyro}}\end{array}\right]$ 18

and in a similar fashion

$p\text{(}{x}_t|u_t,x_{t-1})={\epsilon }_{e_1{\sigma }^2{ }_{\text{vel},\text{gyro}}\left|{\upsilon }_{Gx,\text{gyro}}\right|+e_2{\sigma }^2{ }_{\psi ,\text{gyro}}\left|{\psi }_{\text{gyro}}\right|}\left({\upsilon }_{Gx,\text{actual}}-{\upsilon }_{Gx,\text{gyro}}\right)\cdot {\epsilon }_{e_3{\sigma }^2{ }_{\text{vel},\text{gyro} }\left|{\upsilon }_{Gx,\text{gyro}}\right|+e_4{\sigma }^2{ }_{\psi ,\text{gyro}}\left|{\psi }_{\text{gyro}}\right|}\left(\Delta {\psi }_{\text{actual}}-\Delta {\psi }_{\text{gyro}}\right)$ 19

Environment sensing system and measurement model

The effects of the bicycle’s motion model and road behavior on the perception of the environment and the detection of landmarks that can be used for the mapping procedure and the map updating are multiple. Although they inherently depend on the sensor system deployed to measure its surroundings, the main reasoning behind them remains the same. Thus, we present two cases of sensor designs that can provide environment scans on a two-dimensional level, such as with the use of laser range sensors (LRS), radars, sonars, and so on. Even in the case of three-dimensional sensors, such as multi-level LRS, cameras, or RGB-D cameras, a methodology to extract the distance to the environment objects is used and the same landmark detection limitations that we demonstrate exist, so the main principles of this analysis still apply.

Fixed range finding sensor

At first, we examine the case of a simple range finding sensor attached on the bicycle’s main frame, which is the same as the rear wheel frame for standard bicycles, with the sensor fixed at a position vertical to the ground and its scanning plane horizontal to the ground. As described by Thrun et al., ¹ $p\text{(}{z}_t|x_t,m)$ , the probability of detecting an object k of true range $z_t^{k*}$ with the measurement $z_t^k$ from the pose x_t on the map m, is a mixture of four densities. First, p_hit as a Gaussian of variance σ_hit for small measurement errors around the true range, then, p_short as an exponential distribution with parameter λ_short for errors related to dynamic objects appearing in front of our object, then, p_max as a small uniform density for objects not detected, and finally, p_rand as a uniform distribution for random noise. Their combined distribution can be seen in Figure 2(a), where the maximum range defined by the sensor is z_max. These four different distributions are mixed by a weighted average, defined by the parameters z_hit, z_short, z_max, and z_rand with $z_{\text{hit}}+z_{\text{short}}+z_{\text{max}}+z_{\text{rand}}=1$

$p\text{(}{z}_t^k|x_t,m)={\left[\begin{array}{c}{z}_{\text{hit}}\\ z_{\text{short}}\\ z_{\text{max}}\\ z_{\text{rand}}\end{array}\right]}^{\text{T}}\left[\begin{array}{c}{p}_{\text{hit}}\text{(} z_t^k|x_t,m)\\ p_{\text{short}}\text{(} z_t^k|x_t,m)\\ p_{\text{max}}\text{(} z_t^k|x_t,m)\\ p_{\text{rand}}\text{(} z_t^k|x_t,m)\end{array}\right]$ 20

Figure 2.

“Pseudo-density” of a typical mixture distribution $p( z_t^k|x_t,m)$ of (a) a plain beam model for real object range $z_t^{k*}$ and (b) a likelihood field model with respective nearest objects at distances o₁, o₂, and o₃ as can be seen in (c) (as described in Thrun et al. ¹ ).

If the likelihood field model is used instead of the plain beam model, then we have a two-dimensional representation and the effect of $p_{\text{short}}\text{(} z_t^k|x_t,m)$ is disregarded, but the main behavior remains the same, as can be seen in Figure 2(b). By learning the above model’s intrinsic parameters and having precise knowledge of the sensor’s position and orientation inside the vehicle’s frame, we can directly apply it for the SLAM “update” phase on most ground vehicles carrying a fixed range finding sensor. That is because the scanning level does not suffer any displacement or rotation while the vehicle is moving. However, in the case of an autonomous bicycle, even if the sensor is fixed on the rear wheel frame, the scanning level will be affected by the bicycle’s rotation on the coronal plane, that is, its roll angle ϕ and angular velocity $\dot{\varphi }$ .

As can be seen from Figure 3, given the roll angle ϕ, some obstacles that would normally be detected by the range finding sensor will be missed since their height will be momentarily under the height of the scan plane. On the other hand, from the left part of the figure, we can observe that other obstacles will also be missed, even though they are within the maximum range of the sensor, since the sensor can only observe the road surface after a specific range. In that case, compared to the previous measurement model, even if the parameters of a sensor have been defined, the mixture distribution of the sensor needs to be adjusted. We must incorporate the roll angle ϕ of the bicycle frame, where counterclockwise rotation is positive, meaning ϕ > 0 when the bicycle is tilting to the left, the range of the objects $z_t^k$ , and the bearing of the objects ω, $\omega \in \left(-\pi ,\pi \right]$ , where ω > 0 if an object lies on the bicycle’s right side.

Figure 3.

Range finder scan plane’s rotational effect on object detection measurements.

If $\varphi \cdot \omega >0$ , then we have the first situation described above, which can be referred to as “looking up,” so the entire mixture distribution needs to be multiplied with an exponential distribution of the form $1-e^{\alpha \varphi z_t^k+\beta }$ , where α, β are bicycle-specific, since objects near the bicycle will not be easily missed, compared to further ones. That stops at z_max, leaving probability of z_max unaffected, since for all the objects that will be missed the range finder will return z_max. If $\varphi \cdot \omega <0$ , then we have the “looking down” case, where our sensor may register the road surface. To integrate that effect, we need to calculate the distance of the ground, given by $d_{\text{ground}}= h_{\text{sensor},\varphi }/\left|\text{sin}\varphi \right|$ , where $h_{\text{sensor},\varphi }=h_{\text{sensor}}/\cos \varphi$ is the adjusted sensor height after the bicycle’s roll rotation. Afterwards, we need to multiply the mixture distribution with a uniform unitary filtering distribution, that stops at d_ground if $d_{\text{ground}}<z_{\text{max}}$ , with a height of $z_{\text{max}}/d_{\text{ground}}$ at d_ground, since we want the height to be 1 for $d_{\text{ground}}=z_{\text{max}}$ . That is implemented to account for the fact that if the ground distance d_ground is relatively small, many objects will be lost and the probability of the sensor registering d_ground is high. The above behavior can be observed in Figure 4.

Figure 4.

Adjustment of range finding sensors’ measurement probability distribution when (a) “looking up” and (b) “looking down.”

For 3D SLAM, where measurements are used to construct a three-dimensional representation of the environment, we need to correct the measured position of each object according to the roll angle ϕ at each instant of measurement. The correction transformation should be as follows:

$\left[\begin{array}{c}{x}_{z_t^k,\text{corrected}}\\ y_{z_t^k,\text{corrected}}\\ z_{z_t^k,\text{corrected}}\end{array}\right]=\left[\begin{array}{ccc}1& 0& 0\\ 0& \cos \varphi & -\text{sin}\varphi \\ 0& \text{sin}\varphi & \cos \varphi \end{array}\right]\left[\begin{array}{c}{z}_t^k\cos \omega \\ z_t^k\text{sin}\omega \\ h_{\text{sensor}}/\cos \varphi \end{array}\right]$ 21

Rotating range finding sensor

Since most modern applications require the creation of three-dimensional maps, we need to consider the case of a rotating sensor to cover the maximum field of view (FoV) possible. Most vehicles can incorporate sensors executing complex rotational patterns, covering the entire 360° FoV. For bicycles, however, given their inherent instability, the sensor rotation mechanism and pattern should be designed carefully so that it does not interfere with the balancing motion control of a self-stabilizing bicycle. The most appropriate mounting position is deemed as the frontal part of the bicycle frame, in the position of the headlight. Other favorable positions include the back of the bicycle for monitoring the rear area, or on top of the bicycle, in which case, however, the bicycle cannot be ridden by someone, a trait we wish to maintain. Sensor rotations on the sagittal or the transverse plane are deemed as highly disruptive to the bicycle’s stability, as well as more complex rotations. Therefore, we choose a rotation on the coronal plane, that is, having the range finding sensor tilting up and down. In that way, we can cover the entire frontal and lateral environment of the bicycle, with the only drawback of not being able to sense objects in the bicycle’s rear; for SLAM purposes, however, such a FoV is deemed sufficient.

For two-dimensional measurements, if the sensor’s scanning plane at time t is at a pitch angle θ with the horizontal level (θ = 0), we consider the counterclockwise rotation as positive, that is, θ < 0 when the sensor is facing downwards in the front. Then, the ground in front of the bicycle will be at distance $d_{\text{ground}-\text{front}}=h_{\text{sensor},\theta }/\left|\text{sin}\theta \right|, \theta <0$ , where $h_{\text{sensor},\theta }=h_{\text{sensor}}+r\text{sin}\theta$ is the adjusted sensor height after the sensor’s tilt rotation, if h_sensor is at θ = 0 and r is the radius of the sensor’s tilt rotation, as can be seen in Figure 5(a). If θ < 0 and $\omega \in \left(-\frac{\pi }{2},\frac{\pi }{2}\right)$ or θ > 0 and $\omega \in \left(\pi ,-\frac{\pi }{2}\right)\cup \left(\frac{\pi }{2},\pi \right)$ , then we have a case similar to the first situation in the previous section, where we may miss objects of low height. If θ > 0 and $\omega \in \left(\pi ,-\frac{\pi }{2}\right)\cup \left(\frac{\pi }{2},\pi \right)$ or θ < 0 and $\omega \in \left(-\frac{\pi }{2},\frac{\pi }{2}\right)$ , then we have a case similar to the second situation of the previous section, where we may observe the ground closer than obstacles.

Figure 5.

Ground distance for (a) sensor coronal rotation (tilt) and (b) bicycle roll and sensor tilt integration.

If we now integrate the bicycle roll angle ϕ and the sensor pitch angle θ in one model to calculate the distance of the ground around the bicycle, as can be observed in Figure 5(b), then using the angle that defines the inclination of the entire scanning plane in relation to the ground, that is, the angle at which the scanning plane is vertical to the ground horizon ω_ver

${\omega }_{\text{ver}}={\text{tan}}^{-1}\frac{d_{\text{ground}-\text{front}}}{d_{\text{ground}-\text{side}}}={\text{tan}}^{-1}\frac{h_{\text{sensor}}/\text{sin}\theta }{h_{\text{sensor}}/\text{sin}\varphi }={\text{tan}}^{-1}\frac{\text{sin}\varphi }{\text{sin}\theta }$ 22

we can define the two situations as above. If θ > 0 and $\omega \in \left({\omega }_{\text{ver}}-\frac{\pi }{2},{\omega }_{\text{ver}}+\frac{\pi }{2}\right)$ or θ < 0 and $\omega \notin \left({\omega }_{\text{ver}}-\frac{\pi }{2},{\omega }_{\text{ver}}+\frac{\pi }{2}\right)$ , then we have the first situation described above, “looking up,” where the probability is multiplied with the decreasing exponential, as can be observed in the upper part of Figure 4. If θ < 0 and $\omega \in \left({\omega }_{\text{ver}}-\frac{\pi }{2},{\omega }_{\text{ver}}+\frac{\pi }{2}\right)$ or θ > 0 and $\omega \notin \left({\omega }_{\text{ver}}-\frac{\pi }{2},{\omega }_{\text{ver}}+\frac{\pi }{2}\right)$ , then we have the second situation described above, “looking down,” where the ground distance should be taken into account, as can be seen in the lower part of Figure 4. In this case, the sensor’s height can be calculated as $h_{\text{sensor},\varphi ,\theta }=\left(h_{\text{sensor}}+r\text{sin}\theta \right)/\cos \varphi$ , since we consider the two rotations independent. The different cases of adjusted measurement probability distributions can be observed collectively in Table 1.

Table 1.

Adjusted measurement probability distribution cases.

Rotation	Angular range	Case
Roll	$\varphi \cdot \omega <0$	Looking up
	$\varphi \cdot \omega >0$	Looking down
Roll and pitch	θ > 0 and $\omega \in \left({\omega }_{\text{ver}}-\frac{\pi }{2},{\omega }_{\text{ver}}+\frac{\pi }{2}\right)$ or θ < 0 and $\hspace{1em}\omega \notin \left({\omega }_{\text{ver}}-\frac{\pi }{2},{\omega }_{\text{ver}}+\frac{\pi }{2}\right)$	Looking up
	θ < 0 and $\omega \in \left({\omega }_{\text{ver}}-\frac{\pi }{2},{\omega }_{\text{ver}}+\frac{\pi }{2}\right)$ or θ > 0 and $\hspace{1em} \omega \notin \left({\omega }_{ver}-\frac{\pi }{2},{\omega }_{ver}+\frac{\pi }{2}\right)$	Looking down

Thus, we can see that even if we can reliably control the sensor pitch angle θ, we still need to account for the influence of the bicycle roll angle ϕ, in order to define our landmarks’ or the ground’s detection probability. If we want the observed objects’ corrected position for a three-dimensional representation, we have to apply the following transformation, assuming we have the pitch rotation first and the roll rotation afterwards

$\left[\begin{array}{c}{x}_{z_t^k,\text{corrected}}\\ y_{z_t^k,\text{corrected}}\\ z_{z_t^k,\text{corrected}}\end{array}\right]=\left[\begin{array}{ccc}\cos \theta & \text{sin}\theta \text{sin}\varphi & \text{sin}\theta \cos \varphi \\ 0& \cos \varphi & -\text{sin}\varphi \\ -\text{sin}\theta & \cos \theta \text{sin}\varphi & \cos \theta \cos \varphi \end{array}\right]\left[\begin{array}{c}{z}_t^k\cos \omega \\ z_t^k\text{sin}\omega \\ \frac{h_{\text{sensor}}+r\text{sin}\theta }{\cos \varphi }\end{array}\right]$ 23

In this analysis, it is important to include the specific errors that determine the accuracy of our measurements. Apart from the inherent error inserted by the range finding sensor, which can be represented by the detected range variance ${\sigma }^2{ }_{\text{sensor},z_t^k }$ and the detected bearing variance ${\sigma }^2{ }_{\text{sensor},\omega }$ , we need to consider the rotation error introduced by the motor in charge of the our sensor’s coronal rotation, represented by the tilt angular variance ${\sigma }^2{ }_{\text{motor},\theta }$ . In addition, it is necessary to account for the error introduced by the gyroscopic sensor’s measurement of the bicycle’s roll angle, if a direct measurement is made and we do not use the bicycle’s model, which is represented by the roll angular variance ${\sigma }^2{ }_{\text{gyro},\varphi }$ . Finally, we need to integrate the error in the height of or sensor which depends at each instant on both ${\sigma }^2{ }_{\text{motor},\theta }$ and ${\sigma }^2{ }_{\text{gyro},\varphi }$ .

Implications on SLAM application

Given the above modeling and proposed SLAM framework, we see fit to discuss the limitations of various SLAM algorithms, given the special characteristics of bicycles and akin vehicles, which may cause some algorithms that are potent on other types of vehicles to produce poor results. Then, we can attempt to relate the main uncertainty factor, that of ϕ, to the overall performance of our bicycle SLAM.

At first, as seen in the above analysis, two-wheeled vehicles suffer from an inherent increased degree of uncertainty, a constant systematic, but rather unpredictable clutter, even for relatively simple environment types. This is apparent from the uncertainty of the vehicle’s control, with the need to account for all the errors inserted by the control variables of the motion model or the errors in the odometry measurement, as explained in “Motion and odometry models” section. Moreover, the dependence of the correct detection of landmarks on the rapid and accurate registration of the bicycle’s roll rotation ϕ often results in missed detections or false alarms, as mentioned in “Environment sensing system and measurement model” section. Additional uncertainty is introduced by the sensitivity of such lightweight vehicles to road surface abnormalities. This is fundamentally different from four-wheeled vehicles, where we can rather accurately calculate the vehicle trajectory, given the motion of the wheels or where we can rely on the leveled and rather stable measurements.

Another important aspect of the bicycle SLAM problem is the data association that is often used in various SLAM algorithms that use feature extraction, especially if cameras or RGB-D sensors are used. Many recent feature extraction techniques are becoming more robust to differentiations in viewing angle and perspective, offering scale, rotation, and even affine transform invariant image features. However, the inaccuracy and blurriness introduced by rapid shifts of the sensing plane and by road abnormalities may be detrimental to the application of data association, causing many techniques to produce poor results. That is a common phenomenon for vehicles such as bicycles, where their steering direction and frame roll change rather frequently and very suddenly during their cruise, while avoiding possible obstacles along the way. Such an example can be seen in Figure 6, where the viewing angle of the camera may differ vastly (left–right) and the extracted features from landmarks can vary substantially, with some features not existing in both images or some existing only partially, which can result in great uncertainty of associating the landmarks expected to be registered in a specific location. Given four-wheeled vehicles’ stable motion, even at high speeds, these problems do not have the same gravity, making them easier to overcome.

Figure 6.

Two different views of the same environment under positive and negative bicycle roll angles ϕ. Some landmarks are not in both frames and the sudden shift may create blurry images that inhibit reliable data association.

Finally, given the above uncertainties in localization and environment object detection and recognition, loop closure, an important part in mapping, will suffer from great inaccuracy. Examining moments when bicycle-like vehicles pass from the same spatial positions, having a similar sensing FoV, but different bicycle roll angles ϕ, as can be observed in Figure 7, the environment depiction may differ greatly, registering the ground or missing landmarks as explained in “Environment sensing system and measurement model” section, and the attempted loop closure may even produce undesirable results, something that is not observed in loop closing for vehicles with inherent stability. Thus, loop closure techniques for the bicycle SLAM problem require special attention in their application in order to remain reliable. In general, the selection of a proper SLAM algorithm that can face the above difficulties is the challenge of two-wheeled vehicle SLAM that remains to be examined.

Figure 7.

Loop closure under different bicycle roll angles ϕ, (a) negative, (b) zero, and (c) positive. The landmarks detected are greatly different among these cases, making normal loop closure and map matching unreliable.

Regardless of the specific SLAM algorithm, however, based on our proposed SLAM framework, we can determine the effect that the most important inaccuracy-introducing factor, the bicycle roll angle ϕ, has on the overall SLAM performance. Given the mathematical analysis of the previous sections, we can observe that the greater the variations are in ϕ and $\dot{\varphi }$ , the greater the errors introduced in the calculation of the motion and odometry model will be, resulting in a poorer SLAM performance. It means that for cases and routes that require frequent change in direction, either for closing loops more frequently or for avoiding more obstacles along the path, we expect the SLAM performance to be poorer, producing maps of worse quality, than in the case of more linear routes. That is also expected because in those courses, the more abrupt roll shift and the higher roll angular velocity increase the uncertainty to a degree that makes the bicycle SLAM problem inherently cluttered, regardless of the environmental conditions, and make the landmark matching more difficult. In addition, given the analysis regarding landmark detections of “Environment sensing system and measurement model” section, we expect that trajectories that constantly make the bicycle tilt left or right in specific parts of them, for example, if they involve left or right turns, will produce mappings that include many erroneous landmarks or disregard existing landmarks. That can happen because of the constant registration of the ground as a landmark and due to the sensor missing landmarks that lie outside the sensor’s tilted FoV, respectively.

Autonomous bicycle prototype and experimental setting

To demonstrate the potency of our proposed system to solve the SLAM problem for bicycles and to verify our aforementioned hypotheses, we implement SLAM on the autonomous bicycle prototype designed by our laboratory, as seen in Figure 8(a). The main units of the autonomous bicycle comprise an NI cRIO controller for the control of the throttle, steering, and brake motors, and for the process of the bicycle roll angle data from a gyroscopic sensor, and an Intel NUC mini computer running a combination of Ubuntu and ROS for the collection and process of the perception data from a SICK LiDAR and of the serial input from an Arduino that handles the laser-based speedometer signal. The control system and sensor integration is such that the bicycle maintains its lightweight, rideable, and easily maneuverable structure.

Figure 8.

(a) Autonomous bicycle prototype and experiment venue for (b) type 1 and (c) type 2, 3 experimental procedures.

The implemented self-stabilizing control methodology is described in the work of our team. ¹⁵ We make use of the LiDAR range-bearing data to apply two-dimensional passive SLAM and recreate a mapping of the environment where landmarks are represented by a mixture of weighted Gaussians and the bicycle’s trajectory within it, as a proof-of-concept for our SLAM framework, without involving visual input at this stage. The weights of the estimated map feature Gaussians are determined by the certainty of the landmarks’ existence at that location and are updated after each measurement set is processed. Given its success and recognition as one of the most effective algorithms for SLAM of the last 15 years, we choose to apply the multi-hypothesis (MH) FastSLAM, ²¹ whose multiple hypotheses methodology can theoretically account for the trajectory and measurement uncertainties.

We perform a series of experiments at two different types of experimental venues and conditions. The first, regarding type 1 experiment, involves a rather small terrain of $35\times 15\hspace{0.17em}\text{m}$ with 12 immobile people taking the place of our landmarks, and a total of seven people moving randomly within our terrain, creating a rather noisy environment, as can be seen in Figure 8(b). The random people detections constitute approximately 20% of the overall detections. The second venue, regarding experiment types 2 and 3, consists of an 80 × 80 m baseball field of uneven ground with 31 static black cones as our landmarks, and little other noise caused by moving objects within our terrain, as can be seen in Figure 8(c). In type 2, we make the bicycle run in trajectories consisting of smaller loops, in order to achieve many instances of loop closing, while making sure the bicycle’s course does not include many steering shifts or sudden changes in bicycle roll angle ϕ, resulting in a rather smooth course seen in Figure 9(b). Finally, in type 3, we execute three big clockwise elliptic trajectories around the field while implementing frequent changes in direction, in order to avoid the landmarks/obstacles along the way, which results in a trajectory full of curves, as can be seen in Figure 9(c). In both type 2 and type 3 experiments, we finish our route with a straight crossing of the field to examine if the trajectory estimated by the SLAM algorithm is accurate.

Figure 9.

FastSLAM results for (a) the type 1, (b) type 2, and (c) type 3 experimental cases. The ground truth trajectory is marked with green, the actual landmark locations are marked with green circles, the estimated map features are marked with blue ‘X’s and the estimated bicycle trajectory is marked with red, while odometry registration is marked with yellow.

The bicycle is remotely controlled to execute the desired trajectories and no dynamic path planning or obstacle avoidance techniques are used, in order to ensure the accuracy of these trajectories. All trajectories are registered by our onboard GPS and are used as the ground truth, which can be seen in all the experimental result figures. The routes last 1 − 2 min for type 1 experiment and 4 − 6 min for the types 2 and 3 with the bicycle running at an average speed of 2 m/s. In order to visualize our results, we draw an occupancy grid-like map to depict landmark locations, given that the landmarks’ specific structure is not taken into account. For all experiment types, we draw our estimated landmarks on specific occupancy grid cells according to the location of the corresponding Gaussians’ means, while we choose the occupancy grid’s accuracy as 0.01 m, that is, each grid cell is 1 cm × 1 cm and all landmarks are assumed to occupy only one grid cell, without visualizing their variance. In addition, merging of neighboring Gaussians map features within a specific radius is used in order to ensure nearby detections are related to the same landmark. Moreover, estimated landmarks with low probability of existence, that is, with corresponding Gaussians of low weights under a specific log odds threshold, are pruned after the update phase, in order to ensure the overall map feature count does not increase out of control in case of erroneous detections.

In the experiments, we make use of the velocity-orientation registration odometry model, using the velocity obtained by a rear-wheel laser speedometer and the gyroscopic sensor orientation, and the fixed LiDAR measurement model, meaning that the sensor’s tilt angle θ = 0. Given the bicycle roll angle, the measurement attempts to prevent the registration of ground around the bicycle as an obstacle, by calculating each time the ground horizon, as explained above. After the off-line data process, we construct the most likely map and bicycle trajectory based on the maximum a posteriori technique, accepting the trajectory of the most likely particle at each instant and the map of the most likely particle. For all the experiments, the number of particles used for the MH FastSLAM algorithm was N = 100, while the merging radius is 0.05 m and the pruning log odds threshold is −5.0.

Experimental results

After conducting our experimental procedure, we can verify the capability of our system to provide a reliable SLAM solution. As can be seen by Figure 9(a), where the results of MH FastSLAM are presented for the first type of environment, we observe that despite the increased noise of the environment we can realize SLAM and produce a mapping of the environment. The potency of our methodology is more visible in the second type of experimental venue with experiment types 2 and 3, as shown in Figure 9(b) and (c), respectively, where we have a more accurate mapping of our environment and trajectory representation. In all figures, the ground truth trajectory and actual landmarks are marked with green, while the estimated trajectory is marked with red and the estimated map features are marked with blue ‘X’s.

In all cases, the odometry seems to produce quite accurate results, proving the potency of our framework’s odometry model. The measurement model being used is also capable of predicting the existence or not of the map landmarks with our proposed estimation of the ground distance and the compensation for the missed landmarks. However, especially for the case of the increased environment clutter experiment, we observe that the introduced uncertainty in the registration of the bicycle roll angle ϕ may still cause the registration of the ground as a map feature in multiple occasions or the removal of map features in the case of missed landmarks. This produces numerous false landmark registrations that can cause the trajectory to deviate increasingly and also result in inaccurate mappings.

To provide a quantitative analysis of our experimental results, we also include the feature count of the resulting MH FastSLAM map compared to the real landmark count as they evolve over time for the type 2 and type 3 experimental cases, as can be seen in Figure 10(a) and (c), respectively. The real landmark count increases monotonically as new landmarks enter the vehicle’s FoV. Additionally, we provide the optimal sub-pattern assignment (OSPA) metric ²² for the comparison of the real map with the MH FastSLAM map for those cases, as can be seen in Figure 10(b) and (d). Although the calculation of the OSPA metric is rather complicated, the core equation used is as follows

$\text{OSPA}={\displaystyle \sum _{i=1}^N}{w}_i\sqrt{\frac{SpatErro{r}_i+CardErro{r}_i}{L*Cutoff}}$ 24

where i is the counter for the N particles in the particle set, each having a weight of w_i. SpatError_i is the estimated map feature spatial error of the i-th particle that is calculated as

$SpatError={\displaystyle \sum _L}\{\begin{array}{c}Md*w_G,\hspace{0.17em}\text{if}\hspace{0.17em}Md*w_G<Cutoff\\ Cutoff\text{,}\hspace{0.17em}\text{otherwise}\end{array}$ 25

that is, the sum of the Mahalanobis distances between all the map feature locations and their corresponding closest real landmark locations multiplied with the weight of each feature’s Gaussian, if the map features are estimated within a specific vicinity of the actual landmarks. If the Mahalanobis distance is greater than a threshold called Cutoff, then that threshold is used in the error calculation instead. CardError is the estimated map feature cardinality error that is calculated as

$CardError\hspace{0.17em}=\hspace{0.17em}\text{max}(unaccountedMapWeight,unaccountedEstWeight)\hspace{0.17em}*\hspace{0.17em}Cutoff$ 26

that is, the maximum between the unaccounted actual map weight unaccountedMapWeight, which is greater in case fewer features were estimated to belong to the map than the observable landmarks, and the unaccounted estimated weight unaccountedEstWeight, which is greater in case we have more estimated map features than actual landmarks, multiplied by the Cutoff threshold. Finally, we divide OSPA by the total number of actual observable landmarks L and the Cutoff threshold to normalize the result.

Figure 10.

FastSLAM (a and c) map feature count and (b and d) OSPA results for the type 2 and type 3 experimental cases, respectively.

From the graphs, we can see that despite the inclusion of more landmarks in the maps over the course of time, the map feature count follows the increase of the real landmark count and the OSPA metric is gradually decreasing over time, especially after loop closures. The map feature number increases above the true landmark count due to the erroneous landmark detections mainly from the registration of the ground, as explained in previous sections, leading to large differences at times. This issue is rather evident in Figure 9(c) in the middle of the clockwise elliptic trajectories, where a high number of false landmarks is included in the estimated map, due to the ground registration as the bicycle tilts to the right. This difference decreases when loop closures occur, which can be seen more clearly in Figure 10(c), but still remains relatively high, which constitutes one of the major problems of bicycle SLAM.

Regarding the algorithm’s performance, after the collection of the data, the off-line process time was approximately 5 min for the type 1 experiment and approximately 10 min for types 2 and 3, while the SLAM algorithm was executed on an Intel i5 processor @2.1 GHz over ROS. Although we did not reach online performance with this proof-of-concept implementation, the algorithm’s performance can be improved with a better selection of merging and pruning techniques, as well as a more optimized SLAM algorithm running stand-alone outside ROS, in order to bring it close to real-time execution.

Regarding our hypothesis relating the bicycle roll angle and angular velocity rate of change and the overall SLAM performance, we examine the graphs of Figure 10. By comparing the estimated landmark counts, we see that in the case of type 3 experiment, where we have more direction shifts and more frequent and more sudden bicycle roll changes, the map feature count and landmark count difference still increases and remains high even though the elliptic trajectory passes from the same location three times. The same holds for the OSPA error that increases after the completion of each elliptic loop, although it should be constantly decreasing, since we have continuous loop closure in theory. On the contrary, in type 2 experiment, due to the continuous loop closures and the smoother change rate of bicycle roll angle and angular velocity, the map feature number is closer to the true landmark count during most of the trajectory and gets reduced significantly after the big loop closure, after time t = 300 s. In both experiments though, due to a constant turning towards one direction, more specifically to the right, especially for type 3 experiment, we have an increased number of spurious landmarks detections, due to the ground registration, as predicted in our hypothesis.

Conclusion

In this article, we proposed for the first time a comprehensive theoretical framework for the application of SLAM for autonomous bicycles and similar two-wheeled vehicles. We introduced the motion and odometry models that are specifically applicable for bicycles based on the complex kinematics and dynamics. Furthermore, we examined two different sensor designs that are suitable for bicycles and analyzed the modification of the measurement model for these two cases. Additionally, we presented our laboratory’s autonomous bicycle prototype and the systems that allow it to realize SLAM and implemented two types of experiments based on the MH FastSLAM methodology as a proof-of-concept. The results of our experiments validate the potency of our bicycle SLAM framework, making it clear that our approach can be also applied to similar two-wheeled vehicles, such as motorcycles and others. However, several issues were raised about the specific SLAM algorithm that is applied on top of our framework. In addition, we analyzed the relationship between the rate of change of the bicycle frame roll angle and the overall SLAM performance and verified this hypothesis through our experiments.

We hope that the suggested theoretical framework will provide the proper groundwork for more and more detailed research on the area of autonomous two-wheeled vehicles, especially with regards to environment perception and SLAM. In that way, we hope that autonomous bicycles and akin vehicles can reach the same level of progress as that of autonomous cars. Future work should focus on further testing our proposed framework for more two-wheeled vehicle SLAM applications and on researching SLAM algorithms with which the aforementioned problems of “Implications on SLAM application” section can be faced. The high degree of uncertainty could be dealt with SLAM methodologies that can withstand the high level of clutter and that can preserve the accuracy of the resulting environment map, despite the measurement errors. Secondly, in order to face the problem of data association, methods that include visual odometry or tracking landmarks could be utilized. Finally, to avoid erroneous loop closure instances, methodologies that readjust their maps in either probabilistic ways or by choosing to integrate all the detected objects could be studied. The two-wheeled vehicle SLAM problem is still a new and challenging field of research.