0 引 言
与月球、火星等近地天体相比，小天体的引力场极为微弱，逃逸速度较小，这导致机器人在探测过程中极有可能出现附着反弹、姿态失控等问题[7]。多足机器人被认为是在弱引力环境中移动探测的有前途的解决方案[8-9]，因为它们的腿部夹持可以保持机器人与星表的连接，兼具良好的移动和越障能力。Robert等[10]设计了一个受Oedipus启发的顺应性机器人爪，能够支撑不稳定物体并安全抓住不均匀的岩石表面。受到蟑螂和蜘蛛中观察到的微刺机制的启发，Alan等[11]开发了Spinybot系列壁爬机器人。Kellar等[12]开发的RiSE系列壁爬机器人也采用了类似的微刺设计。Parness等[13]将仿生微刺技术扩展到岩石抓取，提出了一种基于分布式向内抓取技术的全面的微刺锚固机制，使用微刺锚作为空间机器人LEMUR IIB的脚趾。此外，David等[14]提出了一种用于小行星探测的多指微刺夹持装置的原型机器人，其原理与Parness等设计的微刺黏附装置相似。然而，对于大多数小天体的脆弱和多孔表面，仅使用过于强力的夹持方法效果并不理想[15]。Ribeiro等[16-17]提出了一种低反应感知运动规划，通过减少运动反应来提高弱引力下多腿移动的安全性。然而，该方法将小天体星表理想化为了水平表面，对于真实情境中小天体的非结构化地形并不完全适用。
小天体探测的另一大挑战是对星表的先验知识不足，随机性和复杂性的地形限制了机器人的自主行走。尽管目前可以通过卫星掠飞[18-20]等手段获取部分小天体的地表信息，但要完全掌握其细致的地形分布仍然困难，因此大部分小天体星表可以被认为是非结构化地形。通过判定机器人足端与小天体星表的接触状态是适应非结构化地形的一大策略。Gerardo等[21]提出了一种基于接触检测的非结构化地形感知技术，通过力学反馈建立概率模型来判定机器人腿部的触地情况。然而，受限于弱引力条件，机器人腿部与星表之间的相互作用力十分微弱，基于足端力学反馈检测接触状态并不现实。在没有外部力感应的情况下，一个已经成为经典的方法是应用基于机器人广义动量的扰动观测器[22]。Jemin等[23]使用隐藏马尔可夫模型的推断来融合动态、微分运动学和运动学进行接触检测。然而，其假设可以获取加速度测量值，这在实践中可能无法获得。Camurri等[24]使用基于近似逆动力学的力作为输入，以确定足部是否足够稳固，以便用于基速估计的接触概率先验，但忽略了加速度效应。
针对上述研究的不足，本文基于自主研发的六足机器人移动平台Hexapod Rover，提出了一种适用于小天体弱引力环境与非水平地形的六足机器人足端轨迹规划方法，以提升机器人对小天体星表探测的自主性和安全性。
1 弱引力下非力学反馈的地形坡度估计
通过判断机器人足端与小天体表面接触状态来适应复杂地形是一种重要策略。然而，受限于弱引力条件下力学传感器失敏和机器人足端夹持机构设计，通过传统的力学反馈来判断机器人各足端的接触状态是不现实的。因此，采取不依赖于力学反馈的接触检测技术是必要的。
1.1 基于卡尔曼滤波的触地检测
利用六足机器人本体传感器（如惯性传感器、电机编码器等）得到的步态规划及足端高度等信息，可以估计机器人足端与地面的接触概率。
从宏观角度来看，在六足机器人的运动控制保持稳定且不发散的前提下，其实际足端触地状态与步态规划之间的偏差是有限的。因此，可以基于步态规划的先验知识来估计每条腿足端与地面之间的接触概率$P_i^\varphi $ 。同理，在缺少可靠的环境感知设备的前提下，假设足地接触时的足端高度是呈高斯分布的，可以根据足端高度估计第$i$ 条腿足端与地面之间的接触概率$P_i^h$ 。
在缺失足端力学反馈的情况下，为了提高触地检测的可靠性和准确性，可以将步态相位估计的触地概率作为预测触地概率的输入，将足端高度估计的触地概率作为预测的测量，通过卡尔曼滤波融合这两种概率，最终得到触地概率。
标准的卡尔曼滤波预测方程可以表示为以下形式：
式中，$ {{\hat {\boldsymbol{x}}}_{k - 1}} $ 表示在$k - 1$ 时刻的状态估计；$ {{\hat {\boldsymbol{x}}}_{k,k - 1}} $ 表示预测下一时刻的状态；$ {{{\boldsymbol{u}}}_k} $ 表示$k$ 时刻的控制输入；$ {{\boldsymbol{A}}} $ 表示状态转移矩阵；$ {{\boldsymbol{B}}} $ 表示控制输入矩阵；$ {{{\boldsymbol{P}}}_{k,k - 1}} $ 表示预测下一时刻的误差协方差；$ {{\boldsymbol{Q}}} $ 表示过程噪声协方差矩阵。
对于六足机器人，唯一的期望接触状态由步态规划提供，因此将步态规划估计的各足端触地概率$P_i^\varphi $ ($ i = 1,2,3, \cdot \cdot \cdot ,6 $ )作为控制输入$ {{{\boldsymbol{u}}}_k} $ ：
因此，状态转移矩阵$ {{\boldsymbol{A}}} $ 和控制矩阵$ {{\boldsymbol{B}}} $ 分别定义为：
标准的卡尔曼滤波校正方程为：
式中，$ {{{\boldsymbol{K}}}_k} $ 表示$k$ 时刻的卡尔曼增益；$ {{\boldsymbol{H}}} $ 表示观测矩阵；$ {{\boldsymbol{R}}} $ 表示观测噪声协方差矩阵；$ {{{\boldsymbol{z}}}_k} $ 表示$k$ 时刻的观测值。
根据上文所述，观测值$ {{{\boldsymbol{z}}}_k} $ 由足端高度估计的接触概率$P_i^h$ ($ i = 1,2,\cdots,6 $ )组成：
因此，观测矩阵$ {{\boldsymbol{H}}} $ 也可以写成单位矩阵的形式：
最终通过阈值检测规则得到估计的接触状态：
式中，$ {\hat c_k} $ 表示估计的接触状态；$ {\hat c_{k - 1}} $ 是上一个周期估计的接触状态；$ P(c) $ 是经过卡尔曼滤波融合的接触概率；$ {T_h} = 0.6 $ 、$ {T_l} = 0.4 $ 分别是上、下限阈值。
1.2 非水平地形坡度估计
将小天体星表视为随机起伏的三维平面，则任意时刻六足机器人所在平面可由下述参数方程表征：
式中，$ x $ 表示平面上任意一点空间坐标的$ X $ 轴分量，$ y $ 轴表示该点空间坐标的$ Y $ 轴分量，$ z(x,y) $ 表示该点对应的$ Z $ 轴分量，即高度。
根据触地检测的结果，可以获得$ t $ 时刻以及$ t - 1 $ 时刻共6个足地接触点，它们满足以下等式：
式中，$ p_i^x $ 、$ p_i^y $ 、$ p_i^z $ 分别表示第$ i $ 条腿（$ i = 1,2,3, \cdot \cdot \cdot ,6 $ ）触地时的三轴分量。
式（12）共包含6个方程和3个未知数，属于超定方程，可采用最小二乘法进行求解。由于六足机器人足端与地面的接触在时间上是离散发生的，因此基于触地检测的地形平面方程参数必然会发生突变。为了降低这种参数不连续对运动稳定性的影响，设计了一个一阶低通滤波器对估计的地形平面参数做平滑处理：
式中，$ {\hat{[·]}}_{k} $ 表示$ t = {t_k} $ 时刻的参数预测值；$ \alpha $ 为低通滤波系数，权衡坡度估计结果的实时性与连续性，取$ \alpha = 0.3 $ 。
2 基于最小化运动反应的足端轨迹规划
2.1 弱引力下足地摩擦约束
在极端弱引力条件下，需要着重考虑机器人足端力的约束，以防止在移动过程中出现足端打滑或脱附星表的情况。
对于水平平面，要保证六足机器人每条支撑腿的足端均能够在微重力条件下保持附着、不脱离，并且移动时不会打滑，应满足以下关系式：
式中，$ f_i^x $ 、$ f_i^y $ 、$ f_i^z $ 分别表示机器人第$ i $ 条腿的足端力三轴分量；$ f_i^{z,\max } $ 表示足端力在竖直方向允许的最大值；$ \mu $ 表示接触面的摩擦系数。
式（14）也可以转换成矩阵形式的不等式约束：
根据触地检测的结果，可以计算得到机器人所处地形平面相对于水平面的夹角：
式中，$ \phi $ 表示地形平面相对于水平面的翻滚角；$ \theta $ 表示地形平面相对于水平面的俯仰角；$ \varphi $ 表示地形平面相对于水平面的偏航角。由此，地形平面相对于水平面的旋转矩阵$ {{{\boldsymbol{R}}}_{{\mathrm{pla}}}} $ 可以表示为：
在非水平的地形上，对机器人足端力的摩擦锥约束进行如下调整：
2.2 期望落足点规划
采取Raibert启发式[25]对六足机器人摆动腿的落足点进行规划，足端落足点的期望位置可以表示为：
式中，$ {}^w{{\boldsymbol{p}}}_i^e $ 表示期望落足点在世界系中的位置；$ {}^w{{{\boldsymbol{p}}}_{{\text{Raibert}}}} $ 即为Raibert启发项，它可以根据六足机器人的状态信息对落足点进行调整，提升机器人在移动过程中的鲁棒性。$ {}^w{{{\boldsymbol{p}}}_{{\text{Raibert}}}} $ 具有以下形式：
式中，$ {T_{{c_\phi }}} $ 表示一个步态周期中腿部接触相所占的时间；$ {}^w{{{\boldsymbol{v}}}_0} $ 表示六足机器人基座的三维空间速度；$ {}^w{{\boldsymbol{v}}}_0^e $ 表示六足机器人基座的期望速度；$ k $ 是启发系数，用以调整速度误差对落足点的影响，根据经验和实验效果，取$ k = 0.1 $ 。
式（20）得到的$ {}^w{{\boldsymbol{p}}}_i^e $ 通常只有其水平分量具有意义，用$ {}^w{{\boldsymbol{p}}}_i^{e,x} $ 和$ {}^w{{\boldsymbol{p}}}_i^{e,y} $ 分别表示$ {}^w{{\boldsymbol{p}}}_i^e $ 的$ X $ 轴分量和$ Y $ 轴分量。将$ {}^w{{\boldsymbol{p}}}_i^{e,x} $ 和$ {}^w{{\boldsymbol{p}}}_i^{e,y} $ 代入式（11），即可确定第$ i $ 条腿的期望落足点高度为：
2.3 基于贝塞尔曲线的足端轨迹生成
机器人任何部分的移动都会产生反作用力，这可能导致机器人对小天体星表附着失效。本节提出了一种生成摆动腿的足端轨迹的方法，通过最小化机器人的动量变化来降低由摆动腿运动而产生的反作用力。
7阶贝塞尔曲线可表征为如下形式：
式中，$ {{\boldsymbol{x}}}(t) $ 表示机器人腿部摆动过程中足端的轨迹；$ {{{\boldsymbol{a}}}_i} $ 表示多项式的第$ i $ 个系数向量，也称为控制点；$ {t_0} $ 和$ {t_f} $ 分别表示摆动运动的初始和最终时间。
曲线的第一个和最后一个控制点$ {{{\boldsymbol{a}}}_0} $ 和$ {{{\boldsymbol{a}}}_7} $ 可以定义轨迹的初始和最终位置。曲线的第二个和倒数第二个控制点$ {{{\boldsymbol{a}}}_1} $ 和$ {{{\boldsymbol{a}}}_6} $ 可以确定轨迹的初始和最终的第一导数值，即速度。同理，$ {{{\boldsymbol{a}}}_2} $ 和$ {{{\boldsymbol{a}}}_5} $ 可以确定轨迹初始和最终的加速度情况。
由图1 所示，使用6个控制点，可以在轨迹的始末位置处规定约束足端速度和加速度为0，以避免摆动腿在接触星表时产生突然的运动，防止附着反弹：
式中，$ {{{\boldsymbol{x}}}_0} $ 表示腿部足端抬起前的初始位置，由上一状态给定；$ {{{\boldsymbol{x}}}_f} $ 表示腿部足端摆动后的触地点，由期望落足点$ {}^w{{\boldsymbol{p}}}_i^e $ 给定。
凭借一个可以轻松定义特定点的位置、速度和加速度的平滑曲线，继续进行最小化机器人运动反作用力的过程。机器人运动时产生的线性和角动量，可以用以下方程表示：
式中，$ {{\boldsymbol{P}}} $ 和$ {{\boldsymbol{L}}} $ 分别表示机器人总体的线性动量和角动量；$ {m_i} $ 和$ {{{\boldsymbol{I}}}_i} $ 分别表示机器人第$ {\boldsymbol{i}} $ 个杆件的质量和转动惯量；$ {{\dot {\boldsymbol{x}}}_i} $ 和$ {\dot \theta _i} $ 分别表示机器人第$ i $ 个杆件的速度和角速度；$ n $ 表示机器人杆件总数。
通过优化问题计算轨迹的两个剩余系数$ {{{\boldsymbol{a}}}_3} $ 和$ {{{\boldsymbol{a}}}_4} $ ，最小化运动反作用力$ {\dot {\boldsymbol{P}}} $ 和$ {\dot {\boldsymbol{L}}} $ 。目标函数还包括一个关于足端轨迹的竖直高度部分，以保证抬腿离地间隙，避免腿部摆动时与星表碰撞。优化的边界约束是基于机器人配置的关节角度位置，分别对应于最大位置$ {\theta _{\max }} $ 和最小位置$ {\theta _{\min }} $ 。此外，整个过程还需满足足端的摩擦锥约束，以保证机器人对星表的附着。式（27）显示了所谓的优化问题及其约束的定义：
目标函数的各组分的物理含义如下。式（28）和式（29）分别显示了最小化最大运动反作用力的函数$ {J_1}({\dot {\boldsymbol{P}}}(t)) $ 和$ {J_2}({\dot {\boldsymbol{L}}}(t)) $ 。选择最大运动反作用力作为优化指标，因为夹持的脱离通常发生在反作用力达到极限时。如果这个极限没有被克服，即可避免脱离。步高函数$ {J_3}({x_z}(t)) $ 由式（30）描述，确保机器人腿部在最大高度处的尽可能接近期望高度$ h $ ，$ {x_z}(t) $ 表示腿部足端的竖直位置。可调系数$ {k_1} $ 、$ {k_2} $ 和$ {k_3} $ 定义了各优化指标的权重：
3 弱引力下六足机器人移动仿真实验
3.1 仿真对象与仿真环境
本文选取自主研发的六足机器人作为仿真对象，其整体结构与腿部结构如图2 所示。该机器人单腿共包含根关节、髋关节和膝关节三个转动关节，以及基节、股节和胫节三个杆件。表1 ~3 展示了机器人腿部的详细参数[29]。
表 1 腿部杆件物理参数Table 1 Physical parameters of leg members |
| 杆件 | 长度/mm | 质量/g |
| 基节 | 77.57 | 32.0 |
| 股节 | 140.05 | 189.0 |
| 胫节 | 210.38 | 85.0 |
表 2 腿部杆件动力学参数Table 2 Dynamic parameters of leg members |
| 杆件 | 质心位置 （本体系）/mm | 惯量张量/（g·mm2） |
| 基节 | ||
| 股节 | ||
| 胫节 |
表 3 腿部关节参数Table 3 Leg joint parameters |
| 关节 | 转轴（本体系） | 转动范围 | 最大速度/（（°）/s） |
| 根关节 | Z | −30°~30° | 20 |
| 髋关节 | X | −90°~90° | 40 |
| 膝关节 | X | −20°~180° | 40 |
在开源机器人仿真系统Gazebo[30]中，对弱引力条件下的机器人移动控制方案进行了运动仿真。Gazebo能够提供高保真度的物理模拟，准确有效地模拟机器人运动行为，实现机器人与环境的高效碰撞检测与接触力计算。其运动模拟符合物理定律，为感知与移动算法的验证提供有力支撑。机器人运动仿真参数如表4 所示。
表 4 机器人运动仿真参数Table 4 Simulation parameters of robot motion |
| 参数名称 | 符号 | 数值 |
| 接触检测上限阈值 | Th | 0.6 |
| 接触检测下限阈值 | T1 | 0.4 |
| 低通滤波系数 | α | 0.3 |
| 启发系数 | k | 0.1 |
| 目标函数可调系数 | k1 | 0.4 |
| k2 | 0.4 | |
| k3 | 0.2 | |
| 期望步高/m | h | 0.03 |
在Gazebo中构建粗糙表面来模拟典型小天体火卫一（Phobos）[31]星表的风化层，环境参数见表5 。由于小天体非结构化地形样本稀少且获取困难，因此可载入包含非结构化地形信息的RGB图像和灰度图像，以构建非结构化地形。
表 5 仿真环境参数Table 5 Simulation environment parameters |
| 参数名称 | 符号 | 数值 |
| 引力加速度 | ||
| 逃逸速度 | ||
| 温度 | 233 | |
| 滑动摩擦系数 | 1.00 | |
| 软约束力混合参数 | 1.00 | |
| 软误差减小参数 | 0.20 | |
| 接触点动态刚度系数 | ||
| 接触点动态阻尼系数 |
3.2 仿真分析与验证
如图3 所示，六足机器人运动仿真架构分为上下两层。上层架构包括坡度估计器和足端轨迹规划器，其中坡度估计器根据机器人当前状态实时估计当前地形坡度，并将结果输出至轨迹优化模块。轨迹规划模块根据地形以及机器人状态，在线优化每条腿的足端轨迹。下层为机器人底层控制框架，其接收上层输出的足端轨迹作为输入，并生成关节控制指令。该控制框架包括逆运动学模块和比例-微分（Proportional Derivative，PD）控制器，以精确跟踪目标关节轨迹。
令机器人在弱引力仿真环境中前进行走。图4 展示了机器人在仿真环境中的移动快照，红色虚线描绘了机器人的移动轨迹，箭头则指示了机器人的移动方向。实验结果表明，六足机器人能够在不打滑或翻滚的情况下，稳定穿越非结构化地形。在实验中，机器人的平均移动速度为$ 0.01{\text{ m}}/{\text{s}} $ ，这一速度取决于机器人对地形坡度的估计以及足端轨迹生成所耗费的时间。
在仿真过程中，记录了机器人基座的姿态信息，并重点关注基座俯仰角与地面真实数据的对比。图5 (a)展示了机器人移动过程中基座俯仰角随时间的变化情况，其中蓝色实线表示地面实况，绿色实线代表基座俯仰角。实验结果表明，机器人能够在地形坡度不超过$ {15^ \circ } $ 的地势起伏中自如移动，这证明机器人在实验过程中成功适应了地形的随机变化。
为了验证所提出的移动方法能否有效降低机器人受到的反作用力，并因此增加机器人在低重力环境中的行走稳定性，分别采取本文提出的贝塞尔曲线和传统的三角曲线作为足端轨迹，并比较两种情况下机器人的运动表现。为便于分析不同足端轨迹对机器人移动特性的影响，将机器人限制在水平表面进行移动。
此外，还验证了机器人在移动过程中各腿的足端力是否满足摩擦约束，以确保机器人在弱引力条件下能够稳定附着于星表。
定义时间段$\left[ {{t_1},{t_2}} \right]$ 内足端力的平均值为足端力宏观统计量$ \bar f $ ：
表 6 足端力宏观统计量Table 6 Macroscopic statistics of foot force |
| 足端位置 | |||
| 右前腿 | − | ||
| 右中腿 | − | − | |
| 右后腿 | − | − | |
| 左前腿 | |||
| 左中腿 | |||
| 左后腿 | − |
分析表6 可知，机器人各腿的足端力均满足式（18）所描述的摩擦锥约束，这表明机器人在弱引力条件下能够确保对星表的附着。具体而言，对于足端力在$X$ 方向上的分力$ f_i^x $ ，前腿足端呈现正值，中腿足端接近于0，后腿足端呈现负值，表现为前腿“后挫”，后腿“前屈”；对于足端力在$Y$ 方向的分力$ f_i^y $ ，两侧足端呈相反数，表现为“内收”。各腿部状态符合对夹持的刻画。值得注意的是，$ f_i^y $ 在数值上大于$ f_i^x $ ，这是由于机器人沿$X$ 方向移动，夹持主要表现为垂直于运动方向。
4 结 语
本文研究了适用于小天体弱引力环境与非水平地形的六足机器人足端轨迹规划方法。为提高机器人在弱引力条件下的移动稳定性，本文提出了一种低反应足端轨迹规划方法。该方法基于卡尔曼滤波预测机器人的足地接触情况，进而实时估计当前地形坡度；通过分析机器人足端摩擦约束，规划期望落足点，并利用7阶贝塞尔曲线在线优化足端轨迹。通过最小化运动反应，抑制移动过程中的冲击效应，从而降低机器人脱离星表的风险。最后，通过仿真案例验证了该方法在实现六足机器人在弱引力下非结构地形中移动的有效性。