Next Article in Journal
A Customized Extended Kalman Filter for Removing the Impact of the Magnetometer’s Measurements on Inclination Determination
Next Article in Special Issue
Detection and Control Framework for Unpiloted Ground Support Equipment within the Aircraft Stand
Previous Article in Journal
Integrated Controller Design and Application for CNC Machine Tool Servo Systems Based on Model Reference Adaptive Control and Adaptive Sliding Mode Control
Previous Article in Special Issue
Research on Attitude Detection and Flight Experiment of Coaxial Twin-Rotor UAV
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Vibration-Based Recognition of Wheel–Terrain Interaction for Terramechanics Model Selection and Terrain Parameter Identification for Lugged-Wheel Planetary Rovers

1
State Key Laboratory of Robotics and System, Harbin Institute of Technology, Harbin 150001, China
2
State Key Laboratory of Robotics, Shenyang Institute of Automation, Chinese Academy of Sciences, Shenyang 110016, China
3
Institutes for Robotics and Intelligent Manufacturing, Chinese Academy of Sciences, Shenyang 110169, China
4
Key Laboratory of Marine Robotics, Liaoning Province, Shenyang 110169, China
*
Author to whom correspondence should be addressed.
Sensors 2023, 23(24), 9752; https://doi.org/10.3390/s23249752
Submission received: 24 November 2023 / Revised: 8 December 2023 / Accepted: 9 December 2023 / Published: 11 December 2023

Abstract

:
Identifying terrain parameters is important for high-fidelity simulation and high-performance control of planetary rovers. The wheel–terrain interaction classes (WTICs) are usually different for rovers traversing various types of terrain. Every terramechanics model corresponds to its wheel–terrain interaction class (WTIC). Therefore, for terrain parameter identification of the terramechanics model when rovers traverse various terrains, terramechanics model switching corresponding to the WTIC needs to be solved. This paper proposes a speed-independent vibration-based method for WTIC recognition to switch the terramechanics model and then identify its terrain parameters. In order to switch terramechanics models, wheel–terrain interactions are divided into three classes. Three vibration models of wheels under three WTICs have been built and analyzed. Vibration features in the models are extracted and non-dimensionalized to be independent of wheel speed. A vibration-feature-based recognition method of the WTIC is proposed. Then, the terrain parameters of the terramechanics model corresponding to the recognized WTIC are identified. Experiment results obtained using a Planetary Rover Prototype show that the identification method of terrain parameters is effective for rovers traversing various terrains. The relative errors of estimated wheel–terrain interaction force with identified terrain parameters are less than 16%, 12%, and 9% for rovers traversing hard, gravel, and sandy terrain, respectively.

1. Introduction

Planetary rovers need to traverse rough terrains during planetary exploration missions. Wheel–terrain interaction has substantial influence on the mobility of rovers [1,2]. The typical surfaces of Mars consist of hard, sandy, and gravel surfaces [3]. Rovers can move quickly on hard, flat, and non-slippery surfaces. However, they are typically required to traverse over sandy areas, where they suffer from wheel sinkage, or over gravel surfaces where violent vibrations may happen. Hazards originating from the terrain surface itself are called non-geometric hazards. Knowledge of the terrain class and other properties can help improve the mobility of the rover by adjusting its path or control strategy in advance. Recognition of terrain class, terramechanics modeling, and terrain parameter identification are essential to optimize the control of high-mobility rovers. Research progress and existing problems related to terrain classification and parameter identification are discussed in this section.
Terrain classification can be divided into exteroceptive and proprioceptive terrain classification. Exteroceptive terrain classification mainly applies vision or range data [4,5,6,7,8]. For example, Otsu et al. [4] used semi-supervised algorithms instead of supervised algorithms for visual classification systems. Support vector machine (SVM) is used to classify data from different sensors. For the excellent performance of deep learning in object classification and recognition, some studies have utilized deep learning algorithms for terrain classification of planetary rovers, achieving high classification accuracy [6,7]. However, the vision- or laser-based classification methods are usually related to the topmost terrain surface, which may not be the load-bearing surface. An example of this is the Mars terrain, which is covered by a thin layer of drift sand.
Proprioceptive terrain classification has focused on vibration-based terrain classification methods [9,10,11,12,13,14]. Brooks and Iagnemma [9] extracted the power spectral density (PSD) of vibration acceleration to classify the terrain into gravel or sand. In a later work, vibration, vision, and wheel–terrain interaction force were combined to construct a self-supervised terrain classification framework for planetary rovers [10]. However, the terramechanics model used in their paper is not suitable for calculating the interaction force of the wheel and hard/gravel terrain. The terrain parameters used in the terramechanics model are also unidentified. Vicente et al. [11] extracted eleven statistical features of vibration accompanied with the PSD of vibration acceleration to form a new feature vector. An extreme learning machine was used to classify terrains based on the feature vector. Some studies applied deep learning algorithms to classify planetary terrain based on proprioceptive sensors [12,13]. Those methods can reduce classification accuracy caused by inappropriate feature extraction and selection. These methods neglect that the variation of the wheel speed (over a given terrain) changes the amplitude and frequency of the measured wheel vibration signal. If rovers move on the same class of terrain at different speeds, the class of terrain may be incorrectly classified. One solution is to increase the number of trained speeds to obtain vibration features at various speeds [14]. However, there is a drawback that numerous samples at the rover’s multiple moving speeds need to be obtained. Another solution is to extract speed-independent vibration features of the wheel that are proposed in this paper.
Terrain parameters are useful for planetary rovers to maximize the wheel traction or minimize the power consumption [15]. Several methods have been reported to estimate terrain parameters [16,17,18,19]. Iagnemma et al. [16] presented a linearized terramechanics model for estimating the internal cohesion and friction angle of soil. However, many other important terrain parameters such as the cohesive modulus, friction modulus, and sinkage exponent in terramechanics models have not been identified. Hutangkabodee et al. [17] identified the parameters, including the internal friction angle, shearing deformation modulus, and lumped pressure sinkage coefficient, using the Newton–Raphson method based on a simplified terramechanics model. Ding et al. presented an improved terramechanics model with consideration of the slip–sinkage and lug effects. Soil parameters were divided into three independent sets by decoupling a wheel–soil terramechanics model. This model can be used to identify Lunar terrain parameters with the Yutu-2 rover [18]. However, these methods cannot be applicable to rovers traversing firm soil because the assumed terramechanics models do not cover interaction between the rover wheels and hard terrain (HT). The soil of HT does not experience a shear deformation, and there is no shear stress in the wheel–terrain interaction. The soil cohesion coefficient (c) or shearing deformation modulus (K) is not part of the shear characteristic parameters in the shear terramechanics model of the interaction between a wheel and HT. To address this problem, Yang et al. [19] set parameters K and c equal to zero to make a terramechanics model of the interaction between a wheel and soft terrain (ST) suitable for the interaction between the wheel and HT. It is necessary to carry out terrain classification to determine parameters such as K and c in a terramechanics model when calculating interaction force or torque.
Conventional terrain classification classifies terrain classes according to their components to help rovers and mobile robots adjust their control strategy and path planning through which a qualitative judgment of the terrain properties is acceptable. However, an autonomous planetary rover should “know” the terrain properties to calculate the wheel–terrain interaction force to maximize the wheel traction or minimize the power consumption and improve the traversability of the rover during its operational phase. Despite efforts in this regard, the terrain parameter identification of various terrain classes and terramechanics model switching for the interaction between the wheel and soft or hard terrain are still issues to be researched for planetary rovers.
As mentioned above, the following problems are still associated with research works on terrain classification and parameter identification for planetary rovers: (1) a suitable terramechanics model needs to be switched to for rovers moving on various types of terrain. The recognition of wheel–terrain interaction class (WTIC) has been rarely studied. (2) Parameter identification of terramechanics models in various types of terrain (not just sand) needs to be studied. And (3) reported vibration-based speed-independent terrain classification methods require numerous training samples at multiple speeds, which is not feasible for planetary rovers. To solve these problems, the following need to be investigated: (1) how to extract speed-independent vibration features and recognize WTIC to determine the corresponding terramechanics model, especially in the regions where the terramechanics model needs to be switched, and (2) how to identify parameters in the terramechanics model for a planetary rover traversing various types of terrain.
In this paper, the speed-independent vibration features are extracted to recognize WTIC first. Terrain parameters in the terramechanics model that is switched to corresponding to the recognized WTIC are then identified. This paper comprehensively utilizes various sensors installed on the planetary rovers, including accelerometers, cameras, and six-dimensional force/torque sensors, to achieve the recognition of WTIC and identification of terrain parameters. The accelerometer is used to collect wheel vibration acceleration data. Speed-independent vibration features are extracted from the acceleration data to recognize WTIC and select the corresponding terramechanics model. The camera can obtain wheel–terrain interaction images to realize the estimation of wheel sinkage and slip ratio using machine vision technology. The wheel–terrain interaction force information of the rover wheel is acquired using the six-dimensional force/torque sensor. Then, the wheel–terrain interaction force, wheel sinkage and slip ratio are input into the terramechanics model to identify terrain parameters.
The main contributions of this paper are as follows: (1) corresponding to the three wheel–terrain interaction classes (WTICs), this paper proposes three vibration models of lugged wheels and simplifies them to be approximated as periodic functions, showing the difference in wheel vibration frequency and amplitude for different WTICs. Through the dimensionless processing of the amplitude, frequency, mean, and mean square value of the vibration acceleration, four speed-independent vibration features are extracted to reflect the influence of the terrain on the vibration of the wheel. (2) A vibration-based WTIC recognition method without requirement for training data is proposed. Results of single-wheel test experiments show that the simple and effective method has recognized 292 samples of 300 samples correctly, so the recognition accuracy is about 97%. (3) Corresponding to the recognized WTIC, the terramechanics model is switched by determining whether the parameters of the shearing deformation modulus (K) and soil cohesion coefficient © in the model are 0. Terrain parameters in the terramechanics model are identified using the least squares method for rovers traversing various types of terrain. Experiment results obtained using a Planetary Rover Prototype moving on three terrain classes show that the proposed vibration-based WTIC recognition method is effective and has recognized 411 samples of 432 samples correctly, and the terrain parameter identification for various terrains is also effective.
The remainder of this paper is organized as follows. Section 2 discusses the vibration models of the wheel for different WTICs. Based on these models, speed-independent features extracted from the vibration acceleration are described in Section 3. Section 4 presents a method for WTIC recognition by analyzing the range of vibration features for different WTICs. In Section 5, a method is proposed for identifying parameters of the terramechanics model in various types of terrain. The method switches the corresponding terramechanics model according to the recognized WTIC and then uses the model to identify terrain parameters. Section 6 discusses the experiment results of the WTIC recognition and parameter identification for a Planetary Rover Prototype. Finally, some concluding remarks are provided in Section 7.

2. Vertical Vibration for Lugged Wheels of Planetary Rovers

2.1. Vertical Dynamics Model of Wheel Motion

A lugged wheel can be modeled as a mass–spring–damper system shown in Figure 1 [20], which is composed of the mass mw that the wheels of the rover need to load, a spring with a stiffness coefficient kw, and a damper with a damping coefficient cw. Here, y(t) is the vertical output displacement, and q(t) is the base excitation, which is a time-varying vertical displacement of the lowest point of the wheel rim. In addition, Y() and Q() are the Fourier transforms of y(t) and q(t), respectively, where Y() = H()Q(). The frequency response function H() is as follows:
H j ω = c w j ω + k w / m w ω 2 + c w j ω + k w ,
The modulus of H() can be expressed as:
H j ω = 1 + 2 ξ λ 2 / 1 - λ 2 2 + 2 ξ λ 2 ,
where ξ = c w / 2 k w m w is the damping ratio, λ = f v / ω n is the frequency ratio, f v is the frequency of base excitation q(t), and ω n = k w / m w is the natural frequency of the undamped free mass–spring–damper system. According to [20], H j ω is approximately equal to 1 when λ is less than 0.75.
For a lugged wheel on the rover, the spring stiffness kw is very big. The mass mw is assumed as 29 kg, the same as the mass of Spirit rover equally distributed to each wheel. The rover always travels at a low speed, for example, the speed of the Spirit rover is less than 50 mm/s. As a result, the natural frequency ωn is much bigger than the frequency f v of base excitation that is an approximate periodic function for a lugged wheel, as shown later in this paper. Therefore, λ is less than 0.75 and H j ω is approximately equal to 1. Hence,
Y j ω Q j ω ,
According to the properties of a Fourier transform, the following approximation can be obtained.
y t q t .
Therefore, the vertical output displacement and base excitation of the rigid lugged wheel can be regarded as equal in value, i.e., y t = q t . Thus, in the following section, the vibration base excitations of the wheel moving on various terrain classes are analyzed to obtain the wheel vibration output displacements.

2.2. Vibration Models for Wheels Moving on Different Terrains in Wheels’ Center Plane

Wheels with lugs parallel to the wheel axis are usually used on planetary rovers, such as Spirit rover’s wheels. This type of wheel is studied in this paper. The vibrations of this type of wheel result from lug motion and terrain roughness/hardness. Terrains are mainly divided into soft and hard terrain. Both of them can be divided into low-unevenness and high-unevenness terrain. When the planetary rover travels on different classes of terrain, the terramechanics model is different and needs to be switched by setting the values of terrain parameters K and c in the model. In order to determine the values of K and c, wheel–terrain interactions are divided into three classes as shown in Figure 2: lugs entering into soil (LEIS), all lugs contacting terrain surface (ALCT), some lugs contacting terrain surface (SLCT). When the wheel moves on low-unevenness soft terrain and high-unevenness soft terrain, every lug can enter into the soil. When the wheel moves on low-unevenness hard terrain, every lug can contact the terrain surface. When the wheel moves on high-unevenness hard terrain, the same lugs can contact the terrain surface while other lugs cannot contact with the terrain for a wheel rotating 360°. The terrain is classified into three groups: the first class of terrain (T1), which corresponds to WTIC 1: ALCT, prevents each lug from entering the soil during wheel movement, although each lug can be in contact with the terrain. The second class of terrain (T2), which corresponds to WTIC 2: SLCT, makes it impossible for the lugs to enter the soil during the wheel movement, and some lugs cannot be in contact with the terrain. The third class of terrain (T3), which corresponds to WTIC 3: LEIS, refers to the terrain when the wheels transverse it, and where the lugs can enter the soil. The correspondence between terrain classes and the WTICs is shown in Table 1.
The origin of coordinate system Σo-xy is located at the wheel center when the wheel is stationary. The x-axis is parallel to the horizontal direction, and the y-axis is opposite the gravity direction. Σo-xy does not move with the wheel. The angle velocity of the wheel is ωw, l is the lug height, nl is the number of lugs, and αnl (αnl = 2π/nl) is the angle between two adjacent lugs.
According to the number of lugs on the wheel, the wheel can be divided into two classes: olig-lugged and multi-lugged wheels. An olig-lugged wheel is defined as the lugged wheel whose rim can contact the terrain surface when traversing hard terrain (HT). A multi-lugged wheel, which is often used on rovers, is defined as the wheel equipped with a large number of lugs that keep its rim from making contact with the terrain surface when traversing HT.

2.2.1. Vibration Model of Wheels Moving on Flat Hard Terrain

Flat hard terrain is regarded as the reference terrain. The vertical vibration of a wheel traversing the reference terrain is defined as the vertical reference vibration.
Figure 3 shows the movement phases for a wheel rotating around a supporting lug on a reference terrain. Here, θ (−αnl/2 < θαnl/2) is the angle from the direction of the support lug to the gravity direction at the initial time when the wheel starts to move in the anticlockwise direction, which is taken as the positive direction of the angle. At any time, there is no more than one lug in contact with the terrain when an olig-lugged wheel traverses the reference terrain. Assume θ is equal to αnl/2 at the initial time t1 when the wheel starts to move. The movement of the olig-lugged wheel rotating around the support lug can be divided into three phases. (1) Wheel center rising phase: this phase starts at moment t1, when the support lug starts to be in contact with the terrain, until moment t2 when the support lug is perpendicular to the terrain. The wheel center moves from the lowest point to the highest point along the y-axis. (2) Wheel center dropping phase: this phase starts at moment t2 until moment t3 when the wheel rim starts to be in contact with the terrain. During this phase, the wheel center moves from the highest point to the lowest point along the y-axis. (3) Wheel center height unchanging phase: this phase is from moment t3 to moment t4 when the next adjacent lug starts to be in contact with the terrain. The height of the wheel center maintains the lowest invariance.
The wheel repeats the above movement phases with the lugs contacting the terrain in turn. The distance hrt(t) between the lowest point of the wheel rim and the reference terrain surface can be expressed as follows:
h r t ( t ) = max 0 , ( r + l ) cos φ r t ( t ) r .
where r is the wheel radius, l is the lug height, and φ r t ( t ) is the angle between the support lug and gravity direction at the moment t. It can be expressed as
φ r t ( t ) = θ ω w t + α n l × INT ω w t + α n l / 2 θ / α n l .
The   wheel   base   excitation   is   q rt ( t ) = h rt ( t ) min ( h rt ( t ) ) = ( t ) .   Based   on   Equation   ( 5 ) ,   the   vertical   output   displacement   of   the   wheel   is y r t ( t ) = q r t ( t ) = max 0 , ( r + l ) cos φ r t ( t ) r .
For a multi-lugged wheel, the movement of the wheel rotating around the supporting lug can be divided into wheel center rising phase and dropping phase. The wheel center height unchanging phase, in which the wheel rim can be in contact with the terrain, does not exist. The distance hrt(t) is as follows:
h r t ( t ) = r + l cos φ r t ( t ) r .
The output displacement of the wheel is:
y r t t = q r t t = r + l cos φ r t ( t ) cos α n l / 2 r + l .
When the multi-lugged wheel traverses the reference terrain, there are no more than two lugs in contact with the terrain, and the wheel rim cannot be in contact with the terrain surface. Thus, ( r + l ) cos φ r t ( t ) r is bigger than zero. The hrt(t) calculated using Equation (8) is equal to hrt(t) calculated using Equation (5). Equation (8) can be regarded as a special case of Equation (5). Whether the lugs on the wheel are olig or multi, the output displacement for a lugged wheel traversing the reference terrain can be expressed uniformly as follows:
y r t ( t ) = q r t ( t ) = h r t ( t ) min ( h r t ( t ) ) = max 0 , ( r + l ) cos φ r t ( t ) r max 0 , cos α n l / 2 r + l r .

2.2.2. Vibration Model of Wheels with All Lugs Contacting Terrain Surface

In an actual situation, the surface of HT is rugged, as shown in Figure 4. The class of wheel–HT interaction is ALCT. The wheel movement on the terrain is similar to the wheel movement on a reference terrain. The characteristics of the wheel movement are irrelative to the state of the wheel at the initial time. Thus, this paper assumes that the lug starts to contact the terrain at the initial moment t1. The terrain coordinate system O-XY is located at the terrain surface. The X-axis is parallel to the horizontal direction, and the Y-axis coincides with the y-axis. The roughness function of the terrain is expressed as g(X). The distance h(t) between the lowest point of the wheel rim and the terrain surface is
h ( t ) = max g S x ( t ) , r + l cos φ t + g X n .
S x ( t ) = X n r + l sin φ t .
φ ( t ) = θ ω w t + n t 1 α n l .
where Sx(t) is the wheel displacement along the x-axis, n is the number of times the lugs contact the terrain from t1 to the moment t, Xn is the position of the supporting lug on the x-axis.
X 1 = r + l cos θ X n + 1 = X n + a 1 r + l + ω w r ( t e n t s n ) .
a 1 = sin φ t s n + sin φ t e n + α n l .
where t s n is the moment when the supporting lug starts to depart from the terrain, and t e n is the moment when the next lug starts to contact the terrain.
During the movement of the wheel, the position of the lug adjacent to the supporting lug on the X-axis is
X n l ( t ) = X n + r + l × η 1 t     t e ( n 1 ) t t s n X n + a 1 r + l + ω w r ( t t s n ) t s n t t e n .
η 1 t = sin φ t cos α n l 1 + cos φ t sin α n l .
The distance between the lug adjacent to the supporting lug and the terrain surface is
h ˜ t = h t h t .
h t = r + l cos φ t + α n l + g X n l t .
where h t is the sum of the distance between the lug adjacent to the supporting lug and wheel center and terrain roughness of position X n l t on the x-axis.
Then, n can be expressed as follows:
n ( t + ) = INT α t / α n l × INT h / h α t ( n + 1 ) α n l INT α t / α n l + INT h / h α t < ( n + 1 ) α n l .
α t = ω w t + α n l / 2 θ .
where t+ is the moment next to t, and n(0) is equal to 1. The vertical displacement of the wheel is
y t = q t = h t min h t .
For a multi-lugged wheel,
X n l ( t ) = X n + r + l × η 1 t .
h ( t ) = r + l cos φ t + g X n .
h r t ( t ) = r + l cos φ r t ( t ) r .
When a multi-lugged wheel traverses HT, since the wheel rim does not contact the terrain, t s n is equal to t e n and r + l cos φ t + g X n is bigger than g S x ( t ) . Thus, the h(t) calculated using Equation (24) is equal to h(t) calculated using Equation (11). Equation (24) can be regarded as a case of Equation (11). Whether the lugs on the wheel are olig or multi, the output displacement for a lugged wheel traversing HT can be expressed uniformly using Equation (22) with h(t) in Equation (11).

2.2.3. Vibration Model of Wheels with Some Lugs Contacting Terrain Surface

In WTIC 2, SLCT, for example, when the wheel moves on gravel terrain (GT), as shown in Figure 5, the lugs will contact gravel during the movement of the wheel. This leads to some of the lugs hanging in the air without contact with the terrain surface or gravel. The angle between two adjacent supporting lugs is no greater than 90°. The distance between the top point of the i-th lug (lui), which will be in contact with the terrain, and the terrain surface is determined as follows:
h ˜ n + i t = h t h n + i t + r .
h n + i t = r + l cos φ t + i α n l + g X ( n + i ) l t .
where X ( n + i ) l t is the position of the i-th lug (lui), which will be in contact with the terrain on the X-axis.
Moreover, the wheel rim may contact the gravel and no lug will be in contact with the terrain some of the time. The minimum distance between the wheel rim and terrain surface is as follows:
h ˜ w t = h t h w t + r .
h w t = max r x w r r 2 x w 2 + g S x t x w .
where S x t is the displacement of the wheel center along the X-axis, xw is the position of a point on the wheel rim relative to the wheel center along the x-axis. X ( n + i ) l t and S x t can be calculated using the following equations.
X ( n + i ) l t = X w + x l i t   h max = h w t = h ( t ) X n + η i t r + l h max h w t   or   h max = h w t h ( t ) .
S x ( t ) = r sin φ t α w t mod α n l + X w h max = h w t = h ( t ) X n r + l sin φ t h max h w   t or   h max = h w t h ( t ) .
h max = max h t , h n + 1 t , h n + 2 t , , h n + m / 4 t , h w t .
η i t = sin φ ( t ) cos i α n l 1 + cos φ ( t ) sin i α n l .
x l i t = r + l sin φ t + i α n l r sin φ t + α w t mod α n l .
X 0 = r + l cos θ X n + 1 = X n l t e n .
α w t = arcsin x w max 0 / r arcsin X n S x t / r + l .
X w = S x t + x w max 0 .
x w max 0 = arg max x w h w t .
where xwmax0 is the value of xw that maximizes h w t at moment t. The number of wheel lugs in contact with the terrain can be expressed as
iff   INT max h n + 1 t , h n + 2 t , , h n + m / 4 t , h w t / h t = 1 max h n + 1 t , h n + 2 t , , h n + m / 4 t , h w t = h n + i t n t + = n t + i iff   INT max h n + 1 t , h n + 2 t , , h n + m / 4 t , h w t / h t = 1 max h n + 1 t , h n + 2 t , , h n + m / 4 t , h w t = h w t n t + = n t + INT α w t / α nl .
where t+ is the moment next to t, h n + i t is the sum of the distance between the i-th lug (lui) and wheel center and terrain roughness of position X ( n + i ) l t on the x-axis, and m is the number of lugs. The distance h(t+) is
h t + = r cos φ t + α w t + mod α n l + g X w - r h max = h w t = h ( t ) r + l cos φ t + + g X n - r h max h w t   or   h max = h w t h ( t ) .
The distance h(0) between the wheel rim and terrain surface at the initial moment t = 0 is equal to max g S x ( 0 ) , r + l cos φ 0 + g X 0 . The vertical displacement of the vibration output is
y t = q t = h t min h t .
When a multi-lugged wheel moves on GT, some lugs will not be in contact with the surface. Thus, the vertical displacement of the multi-lugged wheel can also be expressed as Equation (41).

2.2.4. Vibration Model of Wheels with Lugs Entering the Soil

In WTIC 3: LEIS, the lugs enter into the soil during wheel movement, for example, when a wheel traverses through sandy terrain (ST), as shown in Figure 6. To facilitate an analysis of the vertical output displacement, assume that a lug starts to enter the soil at the initial moment of the wheel movement. The distance h(t) between the lowest point of the wheel rim and the terrain surface can be expressed as follows:
h t = r r + l cos β ω w t y l t r + l cos β .
where y l ( t ) is the position of the top point on the lug starting to enter the soil on the y-axis, β is the angle between the lug starting to enter the soil and gravity direction. The positions of the top points on two adjacent lugs on the y-axis are expressed as y l 1 ( t ) and y l 2 ( t ) , respectively. They have the relationship y l 2 ( t ) = y l 1 ( t α n l / ω w ) when the wheel–terrain interaction condition, including the terrain properties, slip ratio of the wheel, angular velocity of the wheel, load of the wheel, and so on, remains unchanged. Thus, the position of the top points on the frontmost lug in soil on the y-axis at moment t is y f l ( t ) = y l t INT ω w t / α n l × α n l / ω w . It can be used to calculate the distance h(t) at moment t. Therefore,
h t = r r + l cos β ω w t + INT ω w t / α n l × α n l y l t INT ω w t / α n l × α n l / ω w r + l cos β .
According to Equation (43), h(t) is a periodic function with the period αnl/ωw. The value of h(t) (shown in Figure 7a) for an olig-lugged wheel moving on ST varies as follows during the entire process. When the lug starts to contact the terrain at moment t1, it decreases until moment t2 when the lug is parallel to the y-axis and reduces to the minimum. Subsequently, h(t) increases until moment t3 when the lug departs from the soil and reaches the maximum. Then, it remains constant until the next lug contacts the terrain at moment t4. For a multi-lugged wheel, there is more than one lug excavating the soil at the same time. The value of h(t) varies as shown in Figure 7b during the wheel movement. When the lug starts to contact the terrain at moment t1, h(t) decreases until moment t2 when there is one lug in the soil along the y-axis direction and reduces to minimum. Subsequently, h(t) increases until moment t 3 when there are two adjacent lugs symmetrically about the gravity direction and reaches the maximum. Then, it decreases again until the next lug contacts the terrain at moment t4.
The output displacement of the wheel is
y t = q t = h t min h t .

2.3. Simplification of Vibration Outputs for Wheels Traversing Different Terrains

The vibration models built in Equations (22), (41), and (44) are too complex to analyze the vibration characteristics, thus they need to be simplified. For ALCT, it is known from Equation (20) that the amounts of time each lug contacts the terrain surface is not equal and is influenced by the waves of the terrain. The amplitude of y(t) is slightly larger than that for the wheel traversing the reference terrain surface. However, when the wheel rotates 360° on T1 and reference terrain with the same angular velocity, the amount of time and the number of times of the lugs contact the terrain are the same. Therefore, the time of each lug contacting the terrain can be regarded as approximately the same. In addition, y(t) can be regarded as a periodic function whose amplitude is variable. Compared with the output yrt(t) of the reference vibration, y(t) can be expressed as follows:
y t = C M 1 t × max 0 , ( r + l ) cos φ r t ( t ) r C M 1 t × max 0 , cos α n l / 2 r + l r .
where C M 1 ( t ) ( C M 1 1 ) indicates the increase in the amplitude of y(t) caused by the terrain roughness compared to the amplitude of yrt(t).
For SLCT, the amplitude of y(t) is much bigger than yrt(t) for the wheel traversing reference terrain. The number of contact times between wheel and terrain will be less than the number of lugs when the wheel rotates 360°. So, the frequency of y(t) is lower than nlωw/2π. The output displacement y(t) of the wheel can be expressed as follows:
y t = C M 2 t × max 0 , ( r + l ) cos φ ( t ) r C M 2 t × max 0 , cos α n l / 2 r + l r
φ ( t ) = θ ω w t + α n l × INT ω w t + α n l / 2 θ / α n l × C f ( t ) .
where Cf(t) (Cf(t) > 1) indicates the frequency reduction of y(t) incurred by the terrain roughness when compared to the frequency of yrt(t) for the reference terrain. In addition, C M 2 ( t ) ( C M 2 1 ) indicates the increase in the amplitude of y(t) caused by the terrain roughness compared to the amplitude of yrt(t).
For LEIS, every lug can contact the terrain. According to Equations (43) and (44), the frequency of y(t) is nlωw/2π, equal to the frequency of yrt(t). The amplitude of y(t) is smaller than yrt(t) for the wheel traversing reference terrain. Compared with the output yrt(t) of the reference vibration, y(t) can be expressed as
y t = C M 3 ( t ) × max 0 , ( r + l ) cos φ rt ( t ) r C M 3 ( t ) × max 0 , cos α nl / 2 r + l r .
where C M 3 ( t ) ( C M 3 ( t ) < 1 ) indicates the decrease in the amplitude of y(t) incurred from the lug entering the soil when compared to the amplitude of yrt(t).
Equations (45), (46), and (48) show that the amplitude of y(t) is independent of the angular velocity of the wheel. The value of y(t), which can be approximated as a periodic function, can be expressed as follows:
y t = y ˜ t + C y = C y 0 + i = 1 C y i cos i f a t φ i + C y .
where fa is the frequency of output displacement for the vertical vibration, φi is the phase angle of the i-th harmonic, and C y is the constant of y(t) and is equal to −min(h(t)).

3. Analysis and Extraction of Statistical Vibration Features Independent of the Wheel Velocity

It is difficult to measure vertical displacement of a wheel during movement, whereas vertical acceleration of the wheel is easy to obtain. Therefore, the vibration features of the wheels are extracted from the vertical acceleration signal. The vibration acceleration signal is often described by its statistical features such as mean square value, maximum value, kurtosis, and so on. These statistical features are proportional to wheel speed. If a rover moves on low-unevenness terrain at a high speed, the values of statistical features may be close to that for a rover moving on high-unevenness terrain at a low speed. This can cause low-unevenness terrain to be classified as high-unevenness terrain. In order to solve this problem, new speed-independent features are proposed in this paper.

3.1. Extraction of Vibration Features Independent of Wheel Velocity

According to the wheel’s vertical vibration displacement y(t) in Equation (49), its vertical vibration acceleration can be expressed as a y t = y ¨ t = ω w 2 C y y = ω w 2 C y 0 + i = 1 C y i cos i f a t φ i (where C y is a constant component of y). The vertical vibration acceleration of the wheel can be regard as an approximate period signal, for example, the vertical acceleration of the wheel is an approximate period signal as shown in Figure 8 when the wheel moves on ST. The Fourier-series expansion of a periodic acceleration signal is described as follows:
a y t = C 0 + i = 1 C i cos i f a t φ i .
where Ci is the amplitude of the i-th harmonic.
According to the analysis in Section 2, the changing of wheel–terrain interaction class changes the frequency and amplitude of the vertical vibration acceleration signals of a wheel. The vertical reference acceleration is defined as the vertical acceleration of a wheel when it traverses the reference terrain. Therefore, the vertical reference acceleration can be expressed as a y r t t = y ¨ r t t = ω w 2 C y r t 0 + i = 1 C y r t i cos i f a r t t φ i .
Four statistical features are dimensionlessly processed to be independent of wheel speed. They are as follows:
(1) The frequency ratio of the vertical vibration acceleration for a wheel traversing an actual terrain and the reference terrain is described as
κ f = f a / f a r t = 2 π / n L ω w / 2 π / n L ω w / C ¯ f t = C ¯ f t .
where fa is the frequency of the vertical vibration acceleration for wheels traversing actual terrains with angular velocity ωw, and fart is the frequency of the vertical vibration acceleration for wheels traversing reference terrain with the same angular velocity ωw.
The value of κf is equal to the average value of Cf(t). It is independent of the angular velocity, velocity, and slip ratio of the wheel since Cf(t) is independent of those motion parameters. The value of κf is directly proportional to the frequency of the wheel in contact with the terrain. According to Section 2, the frequency of the vibration signal is appropriately equal to the quotient of the number of times the wheel contacts the terrain divided by the cost time when the wheel rotates 360°. According to Equations (45), (46), and (48), κf is equal to 1 for LEIS and ALCT, and κf is smaller than 1 for SLCT.
(2) The ratio of the sum of all harmonic amplitudes in Equation (50) for a wheel traversing an actual terrain and the reference terrain is described as
κ A = A a / A a r t = ω w 2 C y i / ω w 2 C y r t i = C y i / C y r t i .
where Aa is the sum of all harmonics’ amplitudes of the vertical vibration acceleration for wheels traversing actual terrains with angular velocity ωw, and Aart is the sum of all harmonics’ amplitudes of the vertical vibration acceleration for wheels traversing reference terrain with the same angular velocity ωw.
The κA is independent of the angular velocity, velocity, and slip ratio of the wheel because C y i and C y r t i are independent of those motion parameters. The κA is relative to the stiffness and roughness of the terrain. It is directly proportional to the terrain roughness for terrain that does not cause the wheel to sink and inversely proportional to wheel sinkage for terrain that causes the wheel to sink. For ALCT, κA is slightly greater than 1. For SLCT, κA is much greater than 1. For LEIS, κA is less than 1.
(3) The mean value ratio of vertical vibration acceleration for a wheel traversing an actual terrain and the reference terrain is described as follows:
κ m v = v ¯ a / v ¯ a r t = ω w 2 C y 0 / ω w 2 C y r t 0 = C y 0 / C y r t 0 .
where v ¯ a is the mean value of the vertical vibration acceleration for wheels traversing actual terrains with angular velocity ωw, and v ¯ a r t is the mean value of the vertical vibration acceleration for wheels traversing reference terrain with the same angular velocity ωw.
Due to C y 0 and C y r t 0 being independent of the wheel’s angular velocity, velocity, and slip ratio, the κmv is independent of those motion parameters. It indicates intensity of the vertical speed change.
(4) The mean square ratio of the vertical vibration acceleration for a wheel traversing an actual terrain and the reference terrain is described as follows:
κ m s = s ¯ a / s ¯ a r t = ω w 4 a y 2 t / ω w 4 ¯ / ω w 4 a r t y 2 t / ω w 4 ¯ = a y 2 t / ω w 4 ¯ / a r t y 2 t / ω w 4 ¯ .
where v ¯ a is the mean square value of the vertical vibration acceleration for wheels traversing actual terrains with angular velocity ωw, and v ¯ a r t is the mean square value of the vertical vibration acceleration for wheels traversing reference terrain with the same angular velocity ωw.
The value of κms is independent of the wheel’s angular velocity, velocity, and slip ratio because the numerical value of a y 2 t / ω w 4 is independent of those motion parameters. The κms is relative to the stiffness and roughness of the terrain, similar to κA. It can reflect the variation in speed of the signal from the minimum to the maximum value. It is the powerfulness of vibration energy. For ALCT, κms is slightly greater than 1. For SLCT, κms is much greater than 1. For LEIS, κms is smaller than 1.

3.2. Analysis of Vibration Features Based on a Single-Wheel Experiment

Experiments were conducted using the wheel–terrain interaction test platform with an accelerometer used to collect the vertical acceleration of a wheel, as shown in Figure 9a. The three motors on the platform are the driving, turning, and dragging motors, and the sensors on the platform include the torque, displacement, and six-axis force/torque sensor. The driving motor drives the wheel forward, and the dragging motor is used to simulate the body movement of the rover. A wheel moving with different slip ratios can be controlled by regulating the two motors at different rotational speeds. The displacement sensor is used to measure wheel sinkage. A counterweight is used to adjust vertical load of the wheel, whereas the actual load is measured using the six-axis force/torque sensor. The radius of the wheel is 140 mm and its width is 150 mm. A wheel with 28 lugs of 15 mm in height is shown in Figure 9b. The vertical load of the wheel is 100 N. For each terrain class, the velocities of the wheel are 10, 14, and 20 mm/s, and the slip ratio of the wheel is 0–0.4.
The Mars terrain mainly includes hard, sandy, and gravel surfaces [3]. Rovers need to safely traverse those terrains. Motivated by the challenges of rover mobility, HT, GT, ST, hard terrain with gravel (HGT), and sandy terrain with gravel (SGT) were chosen for the experiment, as indicated in Figure 10. HT is made up of hard slate, GT is constructed with a pile of gravel, ST is sand, HGT is formed by placing gravel randomly on slate, and SGT is sand with a small amount of gravel with a random distribution.
The wheel–HT interaction pertains to ALCT, wheel–GT/HGT interaction pertains to SLCT, and wheel–ST/SGT interaction belongs to LEIS. Figure 11 shows the extraction results for the vibration features. The values of the vibration features are expressed as a fuzzy set {Small (S), Medium (M), Big (B), Very-Big (VB)}. Table 2 shows a comparison of all vibration features for three WTICs. The vibration features extracted to constitute feature vector P for WTIC recognition are {p1, p2, p3, p4} (p1 = lg(κf), p2 = lg(κA), p3 = lg(κmv), and p4 = lg(κms)) According to P, each WTIC can be distinguished from the others.

4. Recognition of Wheel–Terrain Interaction Class Based on Speed-Independent Vibration Features

Every WTIC has a suitable terramechanics model. Based on WTIC recognition, the terramechanics model can be correspondingly selected for identifying terrain parameters. A vibration-based speed-independent WTIC recognition method without requirement of training data is introduced in this section.

4.1. Analysis of Feature Space of the Four Speed-Independent Features

It can be known from Table 2 that the WTIC can be distinguished according to its vibration features’ values. The regions of feature space formed by the vibration features for three WTICs can be partitioned as shown in Figure 12. The value ranges of feature vector P for three WTICs can be known from the analysis of Section 3.2. When the WTIC is ALCT, κf is approximately equal to 1, and κA and κms are slightly bigger than 1. Thus, p1 is approximately equal to 0, and p2 and p4 are slightly bigger than 0. As every lug is in contact with the terrain, there are upper limits for p2 and p4. The experiment results show that p3 is also slightly greater than zero. The feature vector P for ALCT (such as wheel–HT interaction) is located within the vicinity of zero. The vicinity of zero is expressed as {p1, p2, p3, p4|lgεfl < p1 < lgεfu, lgεAl < p2 <lgεAu, lgεmvl < p3 <lgεmvu, lgεmsl < p4 < lgεmsu}, where εfl, εAl, εmvl, and εmsl are the lower limits for κf, κA, κmv, and κms respectively, and εfu, εAu, εmvu, and εmsu are the upper limits for κf, κA, κmv, and κms, respectively.
According to the analysis of Section 3.2, when the WTIC is SLCT, κf is smaller than 1, and κA and κms are greater than 1. The feature vector P for SLCT (such as wheel–HGT/GT interaction) is located in the region {p1, p2, p3, p4 | p1 < lgεfl, lgεAu < p2, lgεmvu < p3, lgεmsu < p4}. When WTIC is LEIS, κf is approximately equal to 1, and κA/κms are smaller than 1. The feature vector P for LEIS (such as wheel–ST/SGT interaction) is in the region {p1, p2, p3, p4 | lgεfl < p1 < lgεfu, p2 < lgεAl, lgεmvl < p3 <lgεmvu, p4 < lgεmsl}. If P is not located in the above three regions, the corresponding WTIC is classified as an unknown wheel–terrain interaction class (UWTIC).

4.2. Analysis of Feature Vector P for the All Lugs Contacting Terrain Surface within the Vicinity of Zero

Owing to the wheels used on the rover all being multi-lugged, the analysis in this section is conducted without considering olig-lugged wheels. When the WTIC is ALCT, the supporting lug of the wheel must depart from the terrain before the adjacent lug along the y-axis. When lug lu1 contacts the terrain after lug lu0 along the y-axis and departs from the terrain before the adjacent lug lu2 along the y-axis, as shown in Figure 13, the vibration displacement caused by lu1 has the highest amplitude in this case. When lu1 is along the y-axis, the distance h(t) between the lowest point of the wheel and the terrain is l. The distance h(t) calculated using lu1 is r + l cos α nl r at the moment when lug lu2 is along the y-axis. The amplitude of the vibration output displacement caused by rotating around lu1 is r + l 1 cos α nl . Thus, to make every lug have contact with the terrain, the amplitude of the vibration output displacement caused by rotating around each lug is no bigger than r + l 1 cos α nl . The amplitude of the reference vibration displacement is r + l 1 cos α nl / 2 . Therefore, C M 1 t in Equation (45) is no bigger than 1 cos α nl / 1 cos α nl / 2 . The following in equalities can be obtained based on Equations (52) and (54).
κ A 1 cos α nl / 1 cos α nl / 2 .
κ m s 1 cos α nl 2 / 1 cos α nl / 2 2 .
Thus, the following equations can be obtained:
ε A u = 1 cos α nl / 1 cos α nl / 2 .
ε m s u = 1 cos α nl 2 / 1 cos α nl / 2 2 .
The parameters εfl, εmvl, εfu, and εmvu are influenced by the noise and accelerometer properties including the sampling frequency and precision. In this study, εfl = 1 − εf, εfu = 1 + εf, εf is set as 0.15, εmvu is set as 10, and the lower limits of εAl, εmvl, and εmsl are 0.5. Therefore, the range of the feature vector P for the ALCT is as follows:
lg 1 ε f p 1 lg 1 + ε f lg ε A l p 2 lg 1 cos α nl / 1 cos α nl / 2 lg ε m v l p 3 lg ε m v u lg ε m s l p 4 lg 1 cos α nl 2 / 1 cos α nl / 2 2 .

4.3. Recognition of Wheel–Terrain Interaction Class Based on Vibration Features

According to the value of feature vector P, the WTIC can be recognized. If P satisfies Equation (59), the WTIC is ALCT. According to the analysis in Section 4.1, the value ranges of feature vector P for SLCT and LEIS are expressed as Equation (60) and Equation (61), respectively.
p 1 < lg 1 ε f lg 1 cos α nl / 1 cos α nl / 2 < p 2 lg ε m v u p 3 lg 1 cos α nl 2 / 1 cos α nl / 2 2 < p 4 .
lg 1 ε f p 1 lg 1 + ε f p 2 < lg ε A l lg ε m v l p 3 lg ε m v u p 4 < lg ε m s l .
If P is not located within the above three regions, the corresponding WTIC is UWTIC. The programming for WTIC recognition is shown in Figure 14. Table 3 shows the WTIC recognition results for the single-wheel test experiments.
The accuracy of the WTIC categorization into the three classes is no less than 96.67%. These results show the ability of the algorithm to distinguish WTICs. When attempting to identify ALCT, the algorithm misclassifies 2 samples of 60 samples as UWTIC, and the recognition accuracy reaches 96.67%. When recognizing SLCT, 2 samples of 120 samples are misclassified as UWTIC, and the recognition accuracy is about 98.34%. Similarly, when recognizing LEIS, 4 samples of 120 samples are misclassified as UWTIC, and the recognition accuracy is approximately 96.67%. The three WTICs are not mutually misclassified. This clearly demonstrates the ability of the algorithm to identify WTICs that induce clearly distinct vibrations on the wheel.

5. Estimation of Terrain Properties Based on the Recognition of Wheel–Terrain Interaction Class

The WTICs can be recognized using the algorithm shown in Figure 14. The terramechanics model is switched corresponding to the recognized WTIC and then used to identify the terrain parameters in this section.

5.1. Switchable Terramechanics Model

Yang et al. [19] inferred that the terramechanics models for wheels moving on sand and hard terrain can be similarly expressed, and the parameters of sand and hard terrain can be unified. The terramechanics model equations for wheel–ST interaction can be switched to the equations that have the same form as the Hunt–Crossley contact model by setting the shearing deformation modulus K and soil cohesion c as zero. It can be used for the interaction between a rigid wheel and hard terrain.
The normal force FN, drawbar pull FDP, and resistance moment MR between the wheel and soil can be expressed as follows:
F N = r b A + B E tan φ k s z N + r b c B E .
F D P = r b A E tan φ B k s z N + r b c A E .
M R = r 2 b c C E + r 2 b C E tan φ k s z N .
where r is the wheel radius; b is the wheel width; z is the wheel sinkage, i.e., deformation of the terrain; c is the soil cohesion; φ is the internal friction angle of the soil; N is the sinkage exponent (N = n0 + n1s); ks is expressed as ks = kc/b + kφ (where kc is the cohesive modulus of the soil and kφ is frictional modulus of soil); and A, B, and C are as follows:
A = cos θ m cos θ 2 θ m θ 2 + cos θ m cos θ 1 θ 1 θ m B = sin θ m sin θ 2 θ m θ 2 + sin θ m sin θ 1 θ 1 θ m C = θ 1 θ 2 2 .
where θ1 (θ1 = arccos(1-z/r)) is the entrance angle, θm (θm = (c1 + c2s)θ1) is the angle of the maximum stress, θ2 (θ2 = c3θ1) is the departure angle, and c1, c2, and c3 are the coefficients of wheel–soil interaction angle. The parameter E is expressed as follows:
E = 1 exp ( j ( θ m ) / K ) .
where K represents the shearing deformation modulus, and j is the shearing displacement. When the wheel contacts hard terrain, the terrain has no shearing deformation. Thus, K is equal to zero. In addition, c can be set to zero according to Yang et al. [19]. Equations (62)–(64) are as follows:
F N = r b A + B tan φ k s z N = k z N .
F D P = A tan φ B A + B tan φ F N = μ F N .
M R = r C tan φ A tan φ B F D P = λ c r F D P .
For interaction between the wheel and hard terrain, θ1 tends toward zero, and thus the coefficient of friction μ tends toward tanφ. The parameter λc is the coefficient used to balance the driving torque. When the parameters K and c are appropriately assigned, the terramechanics model can be used for a rover traversing hard terrain. Therefore, the terramechanics model for wheel–ST interaction can be switched to be suitable for wheel–HT interaction.

5.2. Unified Identification Model for Terrain Parameters

Let X = tan φ E , Y = r b k s z N , and Q = r b c E ; the following equation can be obtained:
F N F D P M R = A 0 B B 0 B A A 0 0 r C r C X Y X Y Q .
Therefore,
X Y + Q = M R / r C .
Y = r C F N B M R / r C A .
F D P = A Y + B ( X Y + Q ) = A 2 + B 2 A C M R r B A F N .
W = F N = A Y + B ( X Y + Q ) = A r b k s z N + B C M R r .
T D = M R = r 2 C E b c + F N tan φ / r A / 1 + B E tan φ / A .
where W and TD are the wheel load and driving torque of the motor, respectively. Equations (73)–(75) are a decoupled terramechanics model regarding unknown terrain parameters Ks, n0, n1, c1, c2, c, φ, and K, which can be used for terrain parameter identification.

5.3. Identification of Terrain Parameters Based on Terramechanics Model Switching

Table 4 shows the correspondence between WTICs and setting values of terrain parameters for terramechanics model switching. WTIC recognition is used to judge whether some parameters should be set to zero to achieve terramechanics model switching. If the WTIC is ALCT, c is set to zero, K is set to zero, E is set to 1. This is because there is no shear deformation or backward accumulation of the soil for ALCT. GT is a typical terrain pertaining to T2 (wheel–T2 interaction is SLCT). When the wheel traverses GT, gravel will move slightly backward and accumulate, which is far less than the accumulation degree of sand. Thus, if the WTIC is SLCT, according to [19], c is set to zero.
The sinkage and slip ratio of the wheel are necessary parameters for identification of the terrain parameters in the terramechanics model. Wheel sinkage can be detected using machine vision methods or through a kinematic analysis of the rover suspension [21,22]. When WTIC is not LEIS, sinkage is set as 0.01mm for terrain parameter identification. The wheel slip ratio can be estimated using a machine vision method [23]. The vertical load W, driving torque TD, and drawbar pull FDP can be measured by installing a six-axis force/torque sensor on the wheel [18].
According to Equation (73), FDP is a function of W, TD, θ1, s, and the unknown parameters c1 and c2. Here, c1 and c2 can be identified according to the measured FDP, W, TD, θ1, and s. According to Equation (74), W is a function of TD, θ1, z, s, identified parameters c1 and c2, and unknown parameters Ks, n0, and n1. The three parameters Ks, n0, and n1 can be identified using the measured W, TD, θ1, z, s, and identified c1, c2. According to Equation (75), TD is a function of W, θ1, s, identified parameters c1 and c2, and unknown parameters c, φ, and K. Parameters c, φ, and K can be identified using the measured W, TD, θ1, s, and identified values of c1 and c2.
The least squares method is used to achieve identification for unknown terrain parameters. When identifying c1 and c2, xe denotes the parameter vector {c1, c1} to be identified; xdata denotes the input data vector {s, FN, TD}; ydata denotes the measured FDP; m denotes the length of xdata and ydata; and F denotes the function of Equation (73). When identifying Ks, n0, and n1, xe denotes the parameter vector {Ks, n0, n1} to be identified; xdata denotes the input data vector {s, TD}; ydata denotes the measured FN; and F denotes the function of Equation (74). When identifying c, φ, and K, xe denotes the parameter vector {c, φ, K} to be identified; xdata denotes the input data vector {s, FN}; ydata denotes the measured TD; and F denotes the function of Equation (75). The ultimate goal is to find the vector xe that best fits Equation (76):
min x e 1 2 F ( x e , x d a t a ) y d a t a = 1 2 i = 1 m [ F ( x e , x d a t a i ) y d a t a i ] 2 .
Figure 15 shows the programming of the terrain parameter identification for different types of terrain.

6. Experiments for Wheel–Terrain Interaction Class Recognition and Terrain Parameter Identification Using a Planetary Rover Prototype

The method proposed in this paper is experimentally validated using a six-wheeled Planetary Rover Prototype. The collected vibration acceleration of the wheel is used to recognize the WTIC. WTIC recognition guides the rover to switch the terramechanics model for identifying terrain parameters in the model using the algorithm shown in Figure 15. The results for the WTIC recognition and identification of terrain parameters in the terramechanics model are discussed in this section.

6.1. Introduction of Planetary Rover Prototype and Terrain Environment

The Planetary Rover Prototype is shown in Figure 16. An accelerator mounted on the front left leg of the rover, near the joint with the wheel axle, is used to collect the vertical acceleration of the wheel to study the vibration-based WTIC recognition. A camera is used to monitor the robot’s wheel and detect the wheel sinkage and slip ratio. The camera is mounted on a supporting arm, with three moving degrees of freedom (DOFs) and three rotational DOFs, installed on the robot body. A gyroscope is attached to the camera to detect the rotation angle of the camera. The position of the camera can be acquired with the relation of the arm and robot body. Every wheel is equipped with a six-axis force/torque sensor to detect the drawbar pull force FDP, vertical load W (equal to normal force FN), lateral force FL, driving moment TD, steering moment MS, and overturning moment MO. The radius of the wheel is 155 mm and its width is 150 mm. The wheel has 24 lugs of 10 mm in height.
Three different surfaces were selected for the experiments: GT, HT (tiled with stone bricks), and ST, as shown in Figure 17. HT pertains to terrain T1, GT pertains to terrain T2, and ST pertains to terrain T3. As such, they represent a wide range of hardnesses from soft to hard, and a roughness ranging from flat to rugged, thus better illustrating the overall accuracy of the methods presented.

6.2. Recognition of Wheel–Terrain Interaction Class for Planetary Rover Prototype

On each of the terrain surfaces, the rover drove back and forth thrice at a constant velocity of 0.003, 0.004, 0.005, 0.006, 0.007, 0.008, 0.009, and 0.01 m/s. During the rover movement, 32-bit vibration acceleration samples were collected using the accelerator shown in Figure 14 at a frequency of 128 kHz. The signal is processed by down-sampling and low-pass filtering to extract vibration features. Figure 18 and Figure 19 show the vibration signal samples for three classes of terrain. Once the single-terrain datasets were collected, multi-terrain datasets were collected with the rover driving from ST to HT to GT in a single traversal, both forward and reverse.
A single dataset consists of vibration data from an about 55 s traversal period. The acceleration signals are segmented with 15 s long intervals to recognize WTIC. Table 5 show the results for the WTIC recognition of single-terrain datasets. They demonstrate an excellent performance of the WTIC recognition method proposed in this paper in distinguishing the WTICs. The results show that 133 samples of 144 wheel–HT interaction (which belongs to ALCT) are recognized as ALCT, and the recognition accuracy reaches 92.36%. Ten samples of wheel–GT interaction (which belongs to SLCT) are misclassified as UWTIC, i.e., 93.06% of the wheel–GT interactions are classified as SLCT. All wheel–ST interaction samples are recognized as LEIS, and the accuracy attains 100%. No misidentification occurred among the three WTICs. The results show that the WTIC recognition algorithm has a high level of accuracy and can effectively distinguish the WTICs.
Figure 20 shows information regarding the vibration acceleration and WTIC recognition results for a distinct single traversal of all three terrain classes shown in Figure 15. The rover is traveling from an area of ST to an area of HT and then on to an area of GT. The data are also divided into 15 s interval segments for WTIC recognition. During the entire movement, the Planetary Rover Prototype operates for a total of 315 s, thus, the acceleration is divided into 21 segments for WTIC recognition. The segmentation is shown in the dotted square in Figure 18, and the WTIC recognition results are represented by the mark at the bottom of the block. These results clearly show the effectiveness of the algorithm in recognizing multiple WTICs during a single traversal. This figure also illustrates the most common error mechanism of the algorithm, which occurs when the rover transitions between homogenous terrain regions.

6.3. Terrain Parameter Identification for Planetary Rover Prototype

The wheel sinkage and slip ratio are necessary parameters for terrain parameter identification and can be detected using machine vision methods shown in our earlier studies [22,23].

6.3.1. Introduction of Sinkage Detection for Wheels Moving on Soft Terrain

Ref. [22] detailed the detection method of wheel sinkage, as shown in Figure 21.
It proposes a new definition of wheel sinkage for rough terrain. Wheel sinkage has been defined as a vector using the wheel–soil boundary. A new vision-based method for sinkage detection is proposed in the paper. The saturation of wheel–terrain image is adjusted using a dynamic piecewise non-linear adjustment to enhance the color contrast between the wheel and soil regions. Through a binarizing process, the wheel section in the color-enhanced image is segmented to acquire an image in which the wheel region is white and other regions are black. Thus, the outline of the wheel region can be extracted from the binary image using the edge detection method. The wheel section outline is corrected to represent the wheel rim operating as a circle. The distance between points of the wheel rim and wheel center is equal to the wheel radius. And the distance between points of the wheel–soil boundary and wheel center is less than the wheel radius. The wheel–soil boundary is extracted based on this characteristic. The sinkage is calculated using the wheel–soil boundary based on the sinkage definition.

6.3.2. Introduction of Wheel Slip Ratio Detection

Ref. [23] detailed a detection method for wheel slip ratio shown in Figure 22. The wheel slip ratio was estimated using the same wheel–terrain images as used in the sinkage detection. At different moments, the wheel position in the wheel–terrain image is the same, while the position of the wheel–trace boundary is different. Thus, the wheel speed can be estimated from the displacement of the wheel trace relative to the wheel at different moments. Based on this idea, Lv et al. [23] built an estimation model for the slip ratio, which is expressed as:
s = 1 Δ S arccos 1 Δ u 2 / 2 r 1 2 r
where ΔS is the displacement of the wheel, Δu is the average displacement of the mark points, and r1 is the radius of the points.
As shown in Figure 22a, the boundaries of the wheel trace were extracted from wheel–terrain images acquired at adjacent moments and were matched to estimate the forward displacement ΔS of the wheel. Meanwhile, the mark points (in this paper, the screws on the wheel are regarded as mark points) of the wheel were also extracted from wheel–terrain images acquired at adjacent moments and were matched to estimate the average displacement Δu of the mark points relative to the wheel center. The slip ratio can then be calculated using Equation (77) according to ΔS and Δu.
When the rover is traversing on HT or GT, there is no obvious wheel trace on the terrain, but there are SURFF feature points that do not exist on ST. The wheel trace boundary extraction and matching proposed by Lv et al. [23] for slip ratio estimation are replaced by the terrain SURFF feature point extraction and matching. The forward displacement ΔS of the wheel can be estimated by detecting the SURFF feature points in the wheel–terrain images acquired at two adjacent moments, thereby estimating the wheel slip ratio, as shown in Figure 22b.

6.3.3. Identification Results of Terrain Parameters for Three Classes of Terrain

The slip ratio of the object wheel can be adjusted by varying the speed of the object wheel and maintaining other wheels at a constant speed. For each terrain class, the ratio is controlled at 0, 0.05, 0.1, 0.15, 0.2, 0.3, 0.4, 0.5, and 0.6. The actual slip ratio is detected visually.
Terrain parameters can be divided into three groups, namely, PI = {c1, c2}, PII = {ks, n0, n1}, and PIII = {c, φ, K}, which are directly related to angles of θm and θ2, normal stress, and shearing stress, respectively. The results of terrain parameter identification using the algorithm in Figure 15 are shown in Table 6. Figure 23 and Figure 24 show comparison of the parameters PII = {ks, n0, n1} and PIII = {c, φ, K} for three terrain classes, respectively, to clearly illustrate the difference in terrain parameters between the three terrain classes.
Wheel–terrain interaction forces are predicted for the wheels without considering the rover dynamics or force coordination among the wheels during the rover movement to simplify the calculation of them for easy control. The estimated drawbar pull force FDP, wheel load W, and driving torque TD with identified terrain parameters are compared with the sensor-measured results as shown in Figure 25, Figure 26 and Figure 27. Their estimated relative errors are shown in Table 7. From Figure 25, Figure 26 and Figure 27 and Table 7, the predicted wheel–terrain interaction forces are close to the measured values. Thus, the simplified predicted wheel–terrain force can be used for optimal control (such as to maximize the wheel traction or minimize the power consumption) for rovers.

6.3.4. Analysis of Parameter Identification Results of Hard Terrain

Table 6 shows that the ks of HT is much larger than that of ST because the internal bonding force of HT is much greater than that of ST and has a large rigidity; thus, the terrain does not undergo substantial deformation when a wheel moves on it. There is no sinkage of the wheel. The sinkage exponents n0 and n1 have different characteristics for three classes of terrain. The n0 of HT is close to 0.5, whereas n1 is approximately zero. Thus, the sinkage exponent N is not affected by wheel slip ratio. The normal force model of the wheel–HT interaction is a non-linear elastic model as shown in Equation (67). The tangential bond of HT is extremely large. During the interaction between the wheel and HT, the soil will not be damaged. There is no shear deformation or backward accumulation for the terrain. Therefore, there is no shear stress in the force of the wheel–terrain interaction, and only normal stress occurs. From Figure 25a, the fitting curve of estimated wheel load W using the identified parameters is extremely close to the fitting curve of the measured data. The relative error of estimated wheel load W is no more than 11.68% as shown in Table 7.
The drawbar pull force can be calculated using the friction model as shown in Equation (68). When the wheel slips, the drawbar pull force of the wheel is equal to the friction between the wheel and the terrain. The parameters c and K are both zero, φ represents the external friction angle of the terrain, and tanφ is the friction coefficient between the terrain and wheel.
Because the wheel sinkage is close to zero, the maximum stress angle θm is calculated as equal to zero according to the coefficients of the wheel–soil interaction angle.
Figure 25a shows that the estimated values of FDP and TD differ greatly from the measured values when slip ratio is 0. And the maximum relative error of the estimated drawbar pull force is 65.8%, as shown in Table 7. However, the estimated FDP and TD are close to measured data except when the slip ratio is 0, as shown in Figure 25b. The relative errors of estimated drawbar pull force and driving torque in Figure 25b are both no more than 15.52%. When the slip ratio is zero, there is no relative movement between the wheel and terrain, and the drawbar pull force (friction force) of the wheel is static friction force. At this time, the actual drawbar pull force (friction force) is smaller than the dynamic friction force calculated according to φ (the friction coefficient between wheel and terrain is tan(φ)). Thus, the errors of the estimated drawbar pull force and the driving torque are large when slip ratio is zero. Figure 25b shows a comparison among the estimated values of the drawbar pull force, the normal load, and driving torque and the sensor measurement values when the slip ratio is bigger than zero. The fitting curve of estimated data is close to the fitting curve of measured data. The relative errors of estimated wheel–HT interaction force with identified terrain parameters are less than 16% when slip ratio is not zero.

6.3.5. Analysis of Parameter Identification Results of Gravel Terrain

Table 6 shows that GT’s ks is close to HT’s, GT has a large rigidity, and the terrain does not undergo substantial deformation when a wheel moves on it. GT’s n0 is close to 0.5 and n1 is about 0.054. The value of n1 is relatively small, but the value of ks is very large. The normal force is sensitive to the sinkage exponent. Therefore, the normal force in Figure 26 has obvious variation.
Under the action of wheels, gravel will slightly move backward and accumulate, which is far less than the accumulation degree of sand. There is shear displacement between pieces of gravel and shear stress in the force of the terrain on the wheel. Therefore, K is not 0 in the shear parameter, and φ represents the external friction angle of the terrain.
Since the wheel sinkage is close to 0, the maximum stress angle calculated from coefficients of the wheel–soil interaction angle is approximately equal to 0.
The drawbar pull force, normal load, and driving torque of the wheel are estimated using Equations (62)–(64) with the identified parameters. Figure 26 shows that the fitting curves of estimated wheel interaction force are close to the fitting curves of measured data. The maximum errors of estimated drawbar pull force, wheel load, and driving torque are 11.85%, 6.67%, and 8.46% respectively.

6.3.6. Analysis for Parameter Identification Results of Sandy Terrain

Table 6 shows that ks of ST is much smaller than that of GT/HT. The bonding force between sand grains is small, and thus the stiffness of sand is low. The wheel easily sinks into the soil when the rover traverses over ST. The n0 of ST is about 0.8169 and n1 is about 1.2957. Thus, the sinkage exponent N is greatly affected by the wheel slip ratio. The normal stress model is a variable exponential non-linear elastic model as shown in Equation (62).
During the interaction between the wheel and terrain, the soil will be displaced, causing the compression under the wheel and accumulation behind the wheel. The terrain places a large amount of shear stress on the wheel, which is not symmetrically distributed along the y-direction of the wheel. The shear stress is also an important factor in calculating the supporting force of the terrain applied to the wheel.
Owing to the obvious shear deformation, both c and K are greater than zero and much larger than those of GT. In addition, φ is the internal friction angle of the soil.
From Figure 27, it can be seen that the relative error of the estimated drawbar pull force, wheel load, and driving torque using Equations (62)–(64) with the identified parameters is small, and the fitting curves of the estimated data are extremely close to the fitting curves of the measured data. Table 7 shows that the maximum relative errors of the estimated drawbar pull force, wheel load, and driving torque are 8.11%, 7.64%, and 7.75%, respectively.
According to the above analysis of experiment results, it can be seen that the relative errors of estimated wheel–terrain interaction force with identified terrain parameters are less than 16%, 12%, and 9% for rovers traversing hard, gravel, and sandy terrain, respectively.

7. Discussion and Conclusions

This paper has been focused on speed-independent vibration-based WTIC recognition. The terramechanics model was switched for terrain parameter identification according to the recognized WTIC. Concluding remarks of this study are provided as follows:
(1) Three vibration models of rigid lugged wheels have been built and analyzed for three WTICs. The three vertical vibration displacements are all approximately periodic functions. Therefore, the vibration acceleration of a rigid lugged wheel was approximated as a periodic signal. Its frequency, amplitude, mean value, and mean square value were extracted as vibration features describing the wheel–terrain interaction response characteristics and processed to make them non-dimensional and independent of the wheel speed.
(2) A vibration-based method for WTIC recognition without requirement of training data is proposed. The regions of the vibration feature space for three WTICs are analyzed. The WTICs are determined by judging which feature space region the wheel–terrain interaction is located in. Experiment results of the single-wheel test show the simple and effective method has recognized 292 samples out of 300 samples correctly, so the recognition accuracy is about 97%.
(3) The proposed WTIC recognition method can help rovers to switch the terramechanics model when traversing hard or soft terrain by modulating the parameters c and K appropriately. Terrain parameters in the switched terramechanics model are identified using the least squares method, and model parameters have been successfully identified for various terrain classes, particularly for hard terrain/rock and soft sandy terrains that are similar to the main terrain types of Mars. The experiment results obtained using a Planetary Rover Prototype moving on the three terrain classes have demonstrated that the WTIC recognition, terramechanics model switching, and parameter identification are effective. The relative errors of estimated wheel–terrain interaction force with identified terrain parameters are less than 16%, 12%, and 9% for rovers traversing hard, gravel, and sandy terrain, respectively.
A vibration-based method for WTIC recognition is proposed in this paper. Previous planetary terrain analysis studies usually identify terrain types, such as ST, GT, HT, SGT, HGT, and so on [4,5,6,7,8]. When rovers move on some of these terrain types (such as ST and SGT), the WTIC is the same. It is difficult to give comprehensive and accurate statistics on the correspondence between terrain type and WTIC. Thus, these methods cannot recognize WTIC, select terramechanics models, and identify terrain parameters. In addition, terrain classification may result in misidentification between different types of terrain. It can lead to incorrect selection of terramechanics models. The WTIC recognition method proposed in this paper does not involve misidentification between different WTICs. It is more conducive to the selection of terramechanics models and the identification of terrain parameters.
In terms of terrain parameter identification, a method for identifying the terrain parameters when the rover moves on various terrains is proposed in this paper. Previous studies on terrain parameter identification (such as [16,17,18]) focused on identifying the terrain parameters of ST without identifying the terrain parameters of other types of terrain. They have not been able to switch and select terramechanics models. If the terramechanics model of wheel–ST interaction is used for terrain parameter identification of other types of terrain, it can lead to incorrect terrain parameter identification results. In terms of identifying the terrain parameters of ST, the method proposed in this paper is similar to [18]. The accuracy of identified terrain parameters is basically the same, and the relative errors of estimated wheel–terrain interaction force with identified terrain parameters are also close. Compared to previous studies on terrain parameter identification, the proposed method can achieve the terrain parameter identification of various types of terrain. It has improved the parameter identification ability.
As future work, the identification of terrain parameters in a terramechanics model in real time for planetary rovers traversing various types of terrain will be further studied. The parameters identified using this method will be used to improve the estimation accuracy of wheel–terrain interaction force in the control model for a planetary rover. And they can be used to calculate the wheel-driving torque for evaluating the terrain traversability [24]. In addition, to improve the autonomous mobility of a rover, a terrain–environment perception system of a planetary rover can be built by combining this study with vision-based terrain classification. The system can help a planetary rover construct a knowledge base of the visual terrain features and terramechanics properties and roughly estimate the terrain parameters’ value ranges in the terramechanics model for terrains that the rover will move on.

Author Contributions

Conceptualization, F.L., N.L., H.G. and L.D.; methodology, F.L. and N.L.; validation, F.L. and N.L.; formal analysis, F.L. and N.L.; investigation, F.L., N.L. and L.D.; data curation, H.G., L.D. and Z.D.; writing—original draft preparation, F.L.; writing—review and editing, N.L., H.G., L.D. and Z.D.; visualization, F.L.; supervision, H.G., L.D., Z.D., H.Y. and Z.L.; funding acquisition, H.G., N.L. and L.D. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by the National Natural Science Foundation of China, grant number 51822502, Foundation for Innovative Research Groups of the Natural Science Foundation of China, grant number 51521003, National Natural Science Foundation of China, grant number 52005122, the Fundamental Research Funds for the Central Universities, grant number HIT.BRETIV.201903, and the “111 Project”, grant number B07018.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

Data are contained within the article.

Conflicts of Interest

The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript; or in the decision to publish the results.

References

  1. Ellery, A. Environment-robot interaction–the basis for mobility in planetary micro-rovers. Robot. Auton. Syst. 2005, 51, 29–39. [Google Scholar] [CrossRef]
  2. Banerjee, S.; Harrison, J.; Furlong, P.M.; Pavone, M. Adaptive meta-learning for identification of rover-terrain dynamics. arXiv 2020, arXiv:2009.10191. [Google Scholar]
  3. Ono, M.; Fuchs, T.J.; Steffy, A.; Maimone, M.; Yen, J. Risk-aware planetary rover operation: Autonomous terrain classification and path planning. In Proceedings of the 2015 IEEE Aerospace Conference, Big Sky, MT, USA, 7–14 March 2015. [Google Scholar]
  4. Otsu, K.; Ono, M.; Fuchs, T.J.; Baldwin, I.; Kubota, T. Autonomous terrain classification with co- and self-Training approach. IEEE Robot. Autom. Lett. 2016, 1, 814–819. [Google Scholar] [CrossRef]
  5. Lu, L.; Ordonez, C.; Collins, E.G.J.; Coyle, E.J.; Palejiya, D. Terrain surface classification with a control mode update rule using a 2D laser, stripe-based structured light sensor. Robot. Auton. Syst. 2011, 59, 954–965. [Google Scholar] [CrossRef]
  6. Huang, G.; Yang, L.; Cai, Y.; Zhang, D.P. Terrain classification-based rover traverse planner with kinematic constraints for Mars exploration. Planet. Space Sci. 2021, 209, 105371. [Google Scholar] [CrossRef]
  7. Guo, J.L.; Zhang, X.Y.; Dong, Y.P.; Huang, B. Terrain classification using mars raw images based on deep learning algorithms with application to wheeled planetary rovers. J. Terramechanics 2023, 108, 33–38. [Google Scholar] [CrossRef]
  8. Halatci, I.; Brooks, C.A.; Iagnemma, K. A study of visual and tactile terrain classification and classifier fusion for planetary exploration rovers. Robotica 2008, 26, 767–779. [Google Scholar] [CrossRef]
  9. Brooks, C.A.; Iagnemma, K. Vibration-based terrain classification for planetary exploration rovers. IEEE Trans. Robot. 2005, 21, 1185–1191. [Google Scholar] [CrossRef]
  10. Brooks, C.A.; Iagnemma, K. Self-supervised terrain classification for planetary surface exploration rovers. J. Field Robot. 2012, 29, 445–468. [Google Scholar] [CrossRef]
  11. Vicente, A.; Liu, J.; Yang, G.Z. Surface classification based on vibration on omni-wheel mobile base. In Proceedings of the IEEE/RSJ International Conference on Intelligent Robots and Systems, Hamburg, Germany, 28 September–2 October 2015. [Google Scholar]
  12. Bai, C.C.; Guo, J.F.; Guo, L.L.; Song, J.L. Deep multi-layer perception based terrain classification for planetary exploration rover. Sensors 2019, 19, 3102. [Google Scholar] [CrossRef] [PubMed]
  13. Ugenti, A.; Vulpi, F.; Domínguez, R.; Cordes, F.; Milella, A.; Reina, G. On the role of feature and signal selection for terrain learning in planetary exploration robots. J. Field Robot. 2022, 39, 355–370. [Google Scholar] [CrossRef]
  14. DuPont, E.M.; Moore, C.A.; Roberts, R.G. Terrain classification for mobile robots traveling at various speeds: An eigenspace manifold approach. In Proceedings of the IEEE International Conference on Robotics and Automation, Pasadena, CA, USA, 19–23 May 2008. [Google Scholar]
  15. Iagnemma, K.; Dubowsky, S. Mobile robot rough-terrain control (RTC) for planetary exploration. In Proceedings of the International Design Engineering Technical Conferences and Computers and Information in Engineering Conference, Baltimore, MD, USA, 10–13 September 2000. [Google Scholar]
  16. Iagnemma, K.; Kang, S.; Shibly, H.; Dubowsky, S. Online terrain parameter estimation for wheeled mobile robots with application to planetary rovers. IEEE Trans. Robot. 2004, 20, 921–927. [Google Scholar] [CrossRef]
  17. Hutangkabodee, S.; Zweiri, Y.H.; Seneviratne, L.D.; Althoefer, K. Soil parameter identification for wheel-soil interaction dynamics and traversability prediction. Int. J. Autom. Comput. 2006, 3, 244–251. [Google Scholar] [CrossRef]
  18. Ding, L.; Zhou, R.; Yuan, Y.; Yang, H.; Li, J.; Yu, T.; Liu, C.; Wang, J.; Li, S.; Gao, H.; et al. A 2-year locomotive exploration and scientific investigation of the lunar farside by the Yutu-2 rover. Sci. Robot. 2022, 7, eabj6660. [Google Scholar] [CrossRef] [PubMed]
  19. Yang, H.G.; Ding, L.; Gao, H.B.; Wang, Z.Y.; Lan, Q.N.; Liu, G.J.; Liu, Z.; Li, W.H.; Deng, Z.Q. High-fidelity dynamic modeling and simulation of planetary rovers using single-input-multi-output joints with terrain property map. IEEE Trans. Robot. 2022, 38, 3238–3258. [Google Scholar] [CrossRef]
  20. Zhou, C. Design and Characteristic Simulation of Automotive Shock Absorbers; Machinery Industry Press: Beijing, China, 2012. [Google Scholar]
  21. Brooks, C.A.; Iagnemma, K.; Dubowsky, S. Visual wheel sinkage measurement for planetary rover mobility characterization. Auton. Robot. 2006, 21, 55–64. [Google Scholar] [CrossRef]
  22. Gao, H.B.; Lv, F.T.; Yuan, B.F.; Li, N.; Ding, L.; Li, N.X.; Liu, G.J.; Deng, Z.Q. Sinkage definition and visual detection for planetary rovers wheels on rough terrain based on wheel–soil interaction boundary. Robot. Auton. Syst. 2017, 98, 222–240. [Google Scholar] [CrossRef]
  23. Lv, F.T.; Gao, H.B.; Li, N.; Ding, L.; Deng, Z.Q.; Liu, G.J. Monocular vision-based estimation of wheel slip ratio for planetary rovers in soft terrain. J. Mech. Eng. 2020, 56, 77–85. [Google Scholar]
  24. Hutangkabodee, S.; Zweiri, Y.H.; Seneviratne, L.D.; Althoefer, K. Performance prediction of a wheeled vehicle on unknown terrain using identified soil parameters. In Proceedings of the 2006 IEEE International Conference on Robotics and Automation, Orlando, FL, USA, 15–19 May 2006. [Google Scholar]
Figure 1. Model of a rigid lugged wheel traversing on any terrain.
Figure 1. Model of a rigid lugged wheel traversing on any terrain.
Sensors 23 09752 g001
Figure 2. Wheel–terrain interaction classes (WTICs), with (a) all lugs contacting terrain surface (ALCT), (b) some lugs contacting terrain surface (SLCT), and (c) lugs entering into soil (LEIS).
Figure 2. Wheel–terrain interaction classes (WTICs), with (a) all lugs contacting terrain surface (ALCT), (b) some lugs contacting terrain surface (SLCT), and (c) lugs entering into soil (LEIS).
Sensors 23 09752 g002
Figure 3. Movement phases of a wheel rotating around the supporting lug when the wheel traverses flat hard terrain.
Figure 3. Movement phases of a wheel rotating around the supporting lug when the wheel traverses flat hard terrain.
Sensors 23 09752 g003
Figure 4. Wheel motion states of rotating around the supporting lug on natural HT.
Figure 4. Wheel motion states of rotating around the supporting lug on natural HT.
Sensors 23 09752 g004
Figure 5. Wheel motion states of rotating around the supporting lug on GT.
Figure 5. Wheel motion states of rotating around the supporting lug on GT.
Sensors 23 09752 g005
Figure 6. Wheel motion states when moving on sandy terrain.
Figure 6. Wheel motion states when moving on sandy terrain.
Sensors 23 09752 g006
Figure 7. The distance h(t) between wheel rim’s lowest point and terrain surface, (a) when an olig-lugged wheel moves on ST, (b) when a multi-lugged wheel moves on ST. The black solid line represents the acceleration signal, the red dashed line represents the time point, and the black dashed line represents the coordinate grid.
Figure 7. The distance h(t) between wheel rim’s lowest point and terrain surface, (a) when an olig-lugged wheel moves on ST, (b) when a multi-lugged wheel moves on ST. The black solid line represents the acceleration signal, the red dashed line represents the time point, and the black dashed line represents the coordinate grid.
Sensors 23 09752 g007
Figure 8. Sample of vibration acceleration when a wheel moves on ST.
Figure 8. Sample of vibration acceleration when a wheel moves on ST.
Sensors 23 09752 g008
Figure 9. Experimental equipment, (a) is the wheel–soil interaction test platform and (b) is the experimental wheel.
Figure 9. Experimental equipment, (a) is the wheel–soil interaction test platform and (b) is the experimental wheel.
Sensors 23 09752 g009
Figure 10. Five terrain classes: (a) HT, (b) GT, (c) ST, (d) HGT, and (e) SGT.
Figure 10. Five terrain classes: (a) HT, (b) GT, (c) ST, (d) HGT, and (e) SGT.
Sensors 23 09752 g010
Figure 11. Values of vibration features for different terrain classes. Vibration features are expressed as the fuzzy set {Small (S), Medium (M), Big (B), Very-Big (VB)}.
Figure 11. Values of vibration features for different terrain classes. Vibration features are expressed as the fuzzy set {Small (S), Medium (M), Big (B), Very-Big (VB)}.
Sensors 23 09752 g011
Figure 12. The regions of feature space formed by the vibration features for different WTICs, with (a) the projection in plane p1p2, (b) the projection in plane p1p3, (c) the projection in plane p1p4, (d) the projection in plane p2p3, (e) the projection in plane p2p4, (f) the projection in plane p3p4.
Figure 12. The regions of feature space formed by the vibration features for different WTICs, with (a) the projection in plane p1p2, (b) the projection in plane p1p3, (c) the projection in plane p1p4, (d) the projection in plane p2p3, (e) the projection in plane p2p4, (f) the projection in plane p3p4.
Sensors 23 09752 g012
Figure 13. The wheel–terrain interaction case in which vibration displacement caused by lu1 has the biggest amplitude.
Figure 13. The wheel–terrain interaction case in which vibration displacement caused by lu1 has the biggest amplitude.
Sensors 23 09752 g013
Figure 14. Flow chart for recognition of wheel–terrain interaction class.
Figure 14. Flow chart for recognition of wheel–terrain interaction class.
Sensors 23 09752 g014
Figure 15. Flow chart for identifying parameters of various types of terrain.
Figure 15. Flow chart for identifying parameters of various types of terrain.
Sensors 23 09752 g015
Figure 16. Planetary Rover Prototype.
Figure 16. Planetary Rover Prototype.
Sensors 23 09752 g016
Figure 17. Terrain environment.
Figure 17. Terrain environment.
Sensors 23 09752 g017
Figure 18. Down-sampled acceleration signals when a wheel traverses (a) ST, (b) HT, and (c) GT.
Figure 18. Down-sampled acceleration signals when a wheel traverses (a) ST, (b) HT, and (c) GT.
Sensors 23 09752 g018
Figure 19. Low-pass-filtered acceleration signals when a wheel traverses (a) ST, (b) HT, and (c) GT.
Figure 19. Low-pass-filtered acceleration signals when a wheel traverses (a) ST, (b) HT, and (c) GT.
Sensors 23 09752 g019
Figure 20. Vibration acceleration and terrain classification results for distinct single traversals of all three terrain classes.
Figure 20. Vibration acceleration and terrain classification results for distinct single traversals of all three terrain classes.
Sensors 23 09752 g020
Figure 21. Flow chart of wheel sinkage detection.
Figure 21. Flow chart of wheel sinkage detection.
Sensors 23 09752 g021
Figure 22. Flow charts of wheel slip ratio detection, with (a) for the wheel traversing ST, and (b) for the wheel traversing GT/HT.
Figure 22. Flow charts of wheel slip ratio detection, with (a) for the wheel traversing ST, and (b) for the wheel traversing GT/HT.
Sensors 23 09752 g022
Figure 23. Comparison of parameters ks, n0, and n1 for three terrains.
Figure 23. Comparison of parameters ks, n0, and n1 for three terrains.
Sensors 23 09752 g023
Figure 24. Comparison of parameters c, φ, and K for three terrains.
Figure 24. Comparison of parameters c, φ, and K for three terrains.
Sensors 23 09752 g024
Figure 25. Wheel–terrain interaction force for Planetary Rover Prototype traversing HT, with (a) s = 0–0.6, and (b) s = 0.05–0.6.
Figure 25. Wheel–terrain interaction force for Planetary Rover Prototype traversing HT, with (a) s = 0–0.6, and (b) s = 0.05–0.6.
Sensors 23 09752 g025
Figure 26. Wheel–terrain interaction force for the Planetary Rover Prototype moving on GT.
Figure 26. Wheel–terrain interaction force for the Planetary Rover Prototype moving on GT.
Sensors 23 09752 g026
Figure 27. Wheel–terrain interaction force for the Planetary Rover Prototype moving on ST.
Figure 27. Wheel–terrain interaction force for the Planetary Rover Prototype moving on ST.
Sensors 23 09752 g027
Table 1. Correspondence between terrain classes and the WTICs.
Table 1. Correspondence between terrain classes and the WTICs.
Terrain ClassWTIC
T1all lugs contacting terrain surface (ALCT)
T2some lugs contacting terrain surface (SLCT)
T3lugs entering into soil (LEIS)
Table 2. The comparison of all vibration features for different WTICs.
Table 2. The comparison of all vibration features for different WTICs.
Vibration FeaturesWTIC
ALCTSLCTLEIS
p1 = lg(κf)MSM
p2 = lg(κA)MBS
p3 = lg(κmv)MBVBMB
p4 = lg(κms)MVBS
Feature vector {p1, p2, p3, p4}{M, M, MB, M}{S, B, VB, VB}{M, S, MB, S}
Annotation: MB = S B.
Table 3. WTIC recognition results of 300 samples for the single-wheel test experiments.
Table 3. WTIC recognition results of 300 samples for the single-wheel test experiments.
Actual WTICRecognized WTIC
ALCTSLCTLEISUWTIC
ALCT (wheel–HT interaction)58002
SLCT (wheel–GT/HGT interaction)011802
LEIS (wheel–ST/SGT interaction)001164
Table 4. Correspondence between terrain class and setting values of terrain parameters for terramechanics model switching.
Table 4. Correspondence between terrain class and setting values of terrain parameters for terramechanics model switching.
WTICTerrain Parameters
Soil Cohesion cShearing Deformation Modulus K
ALCT00
SLCT0Without setting value
LEISWithout setting valueWithout setting value
Table 5. WTIC recognition results for Planetary Rover Prototype test experiments.
Table 5. WTIC recognition results for Planetary Rover Prototype test experiments.
Actual WTICRecognized WTIC
ALCTSLCTLEISUWTIC
ALCT (wheel–HT interaction)1330011
SLCT (wheel–GT/HGT interaction)0134010
LEIS (wheel–ST/SGT interaction)001440
Table 6. Results of terrain parameter identification for Planetary Rover Prototype test experiments.
Table 6. Results of terrain parameter identification for Planetary Rover Prototype test experiments.
Terrainc1c2Ks (KPa/mN)n0n1c (KPa)φ (°)K (mm)
T1 (HT)0.806903.9999 × 1050.61040.000033.340
T2 (GT)004.0053 × 1050.56800.054032.170.6
T3 (ST)0.4424−0.63792.4904 × 1030.81691.296308.629.676.6
Table 7. Relative error of data fitting for FDP, W, and TD.
Table 7. Relative error of data fitting for FDP, W, and TD.
Figure CodeFigure 25aFigure 25bFigure 26Figure 27
FDP65.8%15.2%11.85%8.11%
W11.68%11.68%6.67%7.64%
TD84.83%15.52%8.46%7.75%
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Lv, F.; Li, N.; Gao, H.; Ding, L.; Deng, Z.; Yu, H.; Liu, Z. Vibration-Based Recognition of Wheel–Terrain Interaction for Terramechanics Model Selection and Terrain Parameter Identification for Lugged-Wheel Planetary Rovers. Sensors 2023, 23, 9752. https://doi.org/10.3390/s23249752

AMA Style

Lv F, Li N, Gao H, Ding L, Deng Z, Yu H, Liu Z. Vibration-Based Recognition of Wheel–Terrain Interaction for Terramechanics Model Selection and Terrain Parameter Identification for Lugged-Wheel Planetary Rovers. Sensors. 2023; 23(24):9752. https://doi.org/10.3390/s23249752

Chicago/Turabian Style

Lv, Fengtian, Nan Li, Haibo Gao, Liang Ding, Zongquan Deng, Haitao Yu, and Zhen Liu. 2023. "Vibration-Based Recognition of Wheel–Terrain Interaction for Terramechanics Model Selection and Terrain Parameter Identification for Lugged-Wheel Planetary Rovers" Sensors 23, no. 24: 9752. https://doi.org/10.3390/s23249752

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop