好,现在你已经掌握了从关节空间到末端位姿的映射了,差不多就是一个little baby知道自己的胳膊这个样子转动可以到达那个位置啦,但是你该怎么用你的手去握住奶瓶捏(还好有个瓶,小baby哪要考虑审核呀)?这就是逆运动学的问题,也是实际作业中真实的问题。
讲人话就是:对于有n个关节的机械臂,我们已知从关节空间映射到末端位姿的函数T(q),这个通过DH表或者运动旋量的PoE方式可以简单得到,对一个可达的末端位姿X,求解这个未知数为n维向量的高阶非线性方程组:
T
(
q
)
=
X
T(q)=X
T(q)=X,就是逆运动学的问题了。
很难,就算不是非线性方程,线性方程求解还要考虑n是不是等于6,等于6是个方阵还要考虑是不是奇异矩阵,不是6或者不是线性就更不用说了,难,使人听此凋朱颜。
1. 求解方案
- 迭代
小baby是怎么做到的呢?小baby可不会解方程啊,很简单,小baby生来带着一台高性能计算机,人家先差不多给各个关节一个初始值,够不到就微调,说人话就是迭代,甚至是带反馈的迭代,智能算法~
用指数方程搭配牛顿-拉弗森方法可以轻松迭代,当然这对计算机要求很高,除非不要求实时响应。
想了解可以参考 Morden Robotics, B站上也有配套的公开课。 - 运动学解耦
这是一种很多经典机械臂采用的方案,对一个六轴机械臂,如果末端三个旋转轴相交于一点,也就是设计成一种球形手腕,和基座相连的前三个轴则用来确定球形手腕的中心位置,而那三个相交的旋转轴和前三个关节确定的旋转矩阵共同确定末端位姿。
这就把逆运动学分解(解耦)为两个子问题:先从球形手腕中心位置求解前三个关节的坐标,根据前三个关节确定的旋转矩阵和末端位姿求解出球形手腕对应欧拉角的位姿,再对欧拉角矩阵求逆运算解出最后三个关节变量。
下面主要介绍可以解耦的逆运动学方法:
2. 前三个关节的逆向位置求解
我们已知的是末端位姿,也就是位置
o
6
0
o^0_6
o60 ,和姿态
R
6
0
R^0_6
R60 ,最开始我们利用
o
6
0
o^0_6
o60 可以求解出球形手腕的坐标
o
c
0
o^0_c
oc0 因为
o
6
0
o^0_6
o60 是第六个坐标系的原点,也就是最后一个坐标系的原点,而
o
c
0
o^0_c
oc0 是坐标系3、4、5共同的原点,根据DH约定,我们经常规定最后一个Z轴
z
6
z_6
z6 和
z
5
z_5
z5 共线,那两原点之间只需在
z
6
z_6
z6 方向平移
d
6
d_6
d6 就可以了。
z
6
z_6
z6 方向在基坐标系里怎么表示?这就用到旋转矩阵列向量的本质了,其实
R
6
0
R^0_6
R60 的第三列就是
z
6
z_6
z6 在坐标系0里的向量对应的坐标。
所以我们可以求出球形手腕中心坐标:
o
6
0
=
o
c
0
+
d
6
(
R
[
0
0
1
]
)
⇒
o
c
0
=
o
6
0
−
d
6
(
R
[
0
0
1
]
)
=
[
x
c
y
c
z
c
]
\begin{matrix} o_6^0 = o_c^0 + {d_6}\left( {R\left[ {\begin{matrix} 0\\ 0\\ 1 \end{matrix}} \right]} \right)\\ \\ \Rightarrow o_c^0 = o_6^0 - {d_6}\left( {R\left[ {\begin{matrix} 0\\ 0\\ 1 \end{matrix}} \right]} \right) = \left[ {\begin{matrix} {{{\rm{x}}_c}}\\ {{y_c}}\\ {{z_c}} \end{matrix}} \right] \end{matrix}
o60=oc0+d6
R
001
⇒oc0=o60−d6
R
001
=
xcyczc
对不同的运动学配置(也就是机械臂机构的结构RRR、RRP啥的)有不同的解法,但只要是人设计出来的,求解都比较简单:
2.1 一种RRR结构配置求解
如图:
已知球形手腕
z
c
z_c
zc的位置,求
θ
1
,
θ
2
,
θ
3
\theta_1,\theta_2,\theta_3
θ1,θ2,θ3
我们可以用简单的几何方法求解:
观察俯视图:
显然
θ
1
=
A
tan
2
(
x
c
,
y
c
)
θ
1
=
A
tan
2
(
x
c
,
y
c
)
+
π
\begin{array}{l} {\theta _1} = A\tan 2({x_c},{y_c})\\ \\ {\theta _1} = A\tan 2({x_c},{y_c}) + \pi \end{array}
θ1=Atan2(xc,yc)θ1=Atan2(xc,yc)+π
这两个结果差了180度,如果第二、三个关节不变的话球形手腕的位置就变成中心对称的了,所以实际上这两个不同的第一关节角对应的剩余两个关节角也要对称排布。
这里Atan2(x,y)是一个有两个参数的反正切函数:
if x>=0
Atan2(x,y)=arctan
y
x
\frac{y}{x}
xy
else
Atan2(x,y)=arctan
y
x
\frac{y}{x}
xy+
π
\pi
π
其实就是x(也就是
c
o
s
θ
cos\theta
cosθ )大于零的时候(
−
π
2
<
=
θ
<
=
π
2
-\frac{\pi}{2}<=\theta<=\frac{\pi}{2}
−2π<=θ<=2π)Atan2和arctan返回一样的值(
−
π
2
<
=
θ
<
=
π
2
-\frac{\pi}{2}<=\theta<=\frac{\pi}{2}
−2π<=θ<=2π)
x小于零的时候,返回反正切函数值域加上
π
\pi
π之后的值,这样值域就变成了
[
−
π
2
,
3
π
2
]
[-\frac{\pi}{2},\frac{3\pi}{2}]
[−2π,23π]
对一些具有偏置的机构,第一个关节角也会有两个解,被称为左型、右型手臂位形,差不多这个样子:
对于剩下的
θ
2
,
θ
3
\theta_2,\theta3
θ2,θ3:
就到了我们中学送分大题,解三角形了
有以下结果:
r
=
x
c
2
+
y
c
2
s
=
z
c
−
d
1
cos
θ
3
=
r
2
+
s
2
−
a
2
2
−
a
3
2
2
a
2
a
3
=
D
sin
θ
3
=
±
1
−
D
2
θ
3
=
A
tan
2
(
D
,
±
1
−
D
2
)
θ
2
=
A
tan
2
(
r
,
s
)
−
A
tan
(
a
2
+
a
3
c
3
,
a
3
s
3
)
\begin{array}{l} r = \sqrt {x_c^2 + y_c^2} \\ \\ s = {z_c} - {d_1}\\ \\ \cos {\theta _3} = \frac{{{r^2} + {s^2} - a_2^2 - a_3^2}}{{2{a_2}{a_3}}} = D\\ \\ \sin {\theta _3} = \pm \sqrt {1 - {{\rm{D}}^2}} \\ \\ {\theta _3} = A\tan 2(D, \pm \sqrt {1 - {{\rm{D}}^2}} )\\ \\ {\theta _2} = A\tan 2(r,s) - A\tan ({a_2} + {a_3}{c_3},{a_3}{s_3}) \end{array}
r=xc2+yc2s=zc−d1cosθ3=2a2a3r2+s2−a22−a32=Dsinθ3=±1−D2θ3=Atan2(D,±1−D2)θ2=Atan2(r,s)−Atan(a2+a3c3,a3s3)
ps:r、s在第一个图里画出来了,自己参考,求s的时候用了
z
c
−
d
1
z_c-d_1
zc−d1,这里的
d
1
d_1
d1是DH参数,最好可以一眼就看出来他对应哪里的长度,还是讲个熟练度啊。
psps:这里第三个关节角有两个解,对应了肘部向上和肘部向下两个情况,如果再有个偏执算起来就有4个解了。
至此只用几何方法就把第一个求逆子问题解决了,相当简单和高效。
2.2 作业(大同小异,不想写了 : ) )
对于另一种球坐标形式的机械臂,你要不试试怎么求逆(md,体会到老师留作业的快感了)
3. 逆向姿态
回到我们的问题
知道末端位姿,也就是位置和姿态,末端的位置简化为了求球形手腕的位置对应的前三个关节变量现在,我们已经求出来了
θ
1
,
θ
2
,
θ
3
\theta_1,\theta_2,\theta_3
θ1,θ2,θ3, 所以通过DH法我们可以求出
R
3
0
R^0_3
R30
现在的问题就是:已知
R
6
0
R^0_6
R60 和
R
3
0
R^0_3
R30, 求球形手腕的三个欧拉角:
三个欧拉角决定的旋转矩阵其实就是
R
6
3
R^3_6
R63 ,所以:
R
6
0
=
R
3
0
R
6
3
得到:
R
6
3
=
(
R
3
0
)
T
R
6
0
\begin{matrix} R^0_6=R^0_3R^3_6 \\ \\ \text{得到:}R^3_6=(R^0_3)^TR^0_6 \end{matrix}
R60=R30R63得到:R63=(R30)TR60
R
6
3
R^3_6
R63 如果是ZYZ欧拉角生成的话(我见到的好像大部分都是ZYZ),我们把老朋友搬过来:
R
=
[
c
ϕ
c
θ
c
ψ
−
s
ϕ
s
ψ
−
c
ϕ
c
θ
s
ψ
−
s
ϕ
c
ψ
c
ϕ
s
θ
s
ϕ
c
θ
c
ψ
+
c
ϕ
c
ψ
−
s
ϕ
c
θ
s
ψ
+
c
ϕ
c
ψ
s
ϕ
s
θ
−
s
θ
c
ψ
s
θ
s
ψ
c
θ
]
{\rm{R}} = \left[ {\begin{matrix} {{c_\phi }{c_\theta }{c_\psi } - {s_\phi }{s_\psi }}&{ - {c_\phi }{c_\theta }{s_\psi } - {s_\phi }{c_\psi }}&{{c_\phi }{s_\theta }}\\ {{s_\phi }{c_\theta }{c_\psi } + {c_\phi }{c_\psi }}&{ - {s_\phi }{c_\theta }{s_\psi } + {c_\phi }{c_\psi }}&{{s_\phi }{s_\theta }}\\ { - {s_\theta }{c_\psi }}&{{s_\theta }{s_\psi }}&{{c_\theta }} \end{matrix}} \right]
R=
cϕcθcψ−sϕsψsϕcθcψ+cϕcψ−sθcψ−cϕcθsψ−sϕcψ−sϕcθsψ+cϕcψsθsψcϕsθsϕsθcθ
我们接下来的任务就是求三个欧拉角
ϕ
,
θ
,
ψ
\phi,\theta,\psi
ϕ,θ,ψ
下文中R中的元素用
r
i
j
r_{ij}
rij 来表示:
- case1:
当 θ \theta θ不为零的时候
θ = A tan 2 ( r 33 , 1 − r 33 2 ) , 或 θ = A tan 2 ( r 33 , − 1 − r 33 2 ) \theta = A\tan 2({r_{33}},\sqrt {1 - r_{33}^2} ),\text{或} \;\;\theta = A\tan 2({r_{33}}, - \sqrt {1 - r_{33}^2} ) θ=Atan2(r33,1−r332),或θ=Atan2(r33,−1−r332)
-
s
θ
>
0
s_\theta>0
sθ>0
这时:
ϕ = A tan 2 ( r 13 , r 23 ) ψ = A tan 2 ( − r 31 , r 32 ) \begin{array}{l} \phi = A\tan 2({r_{13}},{r_{23}})\\ \\ \psi = A\tan 2( - {r_{31}},{r_{32}}) \end{array} ϕ=Atan2(r13,r23)ψ=Atan2(−r31,r32) -
s
θ
<
0
s_\theta<0
sθ<0
这时
ϕ = A tan 2 ( − r 13 , − r 23 ) ψ = A tan 2 ( r 31 , r 32 ) \begin{array}{l} \phi = A\tan 2({-r_{13}},{-r_{23}})\\ \\ \psi = A\tan 2( r_{31},{r_{32}}) \end{array} ϕ=Atan2(−r13,−r23)ψ=Atan2(r31,r32)
区分这两种情况是因为Atan2的算法中第一个参数的正负会影响结果的区间,所以 r 13 , r 31 r_{13},r_{31} r13,r31的符号会受 s θ s_\theta sθ影响,所以分两个情形进行修正。
- case2
当 θ \theta θ为零的时候
这时候旋转矩阵退化为
[ c ϕ + ψ − s ϕ + ψ 0 s ϕ + ψ c ϕ + ψ 0 0 0 1 ] 或者 [ − c ϕ − ψ − s ϕ − ψ 0 s ϕ − ψ c ϕ − ψ 0 0 0 − 1 ] \left[ {\begin{matrix} {{c_{\phi + \psi }}}&{ - {s_{\phi + \psi }}}&0\\ {{s_{\phi + \psi }}}&{{c_{\phi + \psi }}}&0\\ 0&0&1 \end{matrix}} \right] \text{或者}\ \ \left[ {\begin{matrix} { - {c_{\phi - \psi }}}&{ - {s_{\phi - \psi }}}&0\\ {{s_{\phi - \psi }}}&{{c_{\phi - \psi }}}&0\\ 0&0&{ - 1} \end{matrix}} \right] cϕ+ψsϕ+ψ0−sϕ+ψcϕ+ψ0001 或者 −cϕ−ψsϕ−ψ0−sϕ−ψcϕ−ψ000−1
分别可以解出 ϕ + ψ 和 ϕ − ψ \phi+\psi \text{和} \ \ \phi-\psi ϕ+ψ和 ϕ−ψ的值
这时候有无穷多解
至此逆运动学结束。