吴鹏，洪娟，陈广，徐博

(1.哈尔滨工程大学 自动化学院，黑龙江 哈尔滨 150001； 2.北京航天控制仪器研究所，北京 100854)

串联冗余自由度机构关节角位置的规划算法

吴鹏1，洪娟2，陈广1，徐博1

(1.哈尔滨工程大学 自动化学院，黑龙江 哈尔滨 150001； 2.北京航天控制仪器研究所，北京 100854)

为了建立一种实时求解串联冗余自由度机构关节角位置的方法，解决串联冗余自由度机构逆运动学求解中角位置的求解问题，给出了逆运动学的数值解。根据加权伪逆理论，采用瞬时最优速度表达关节空间角增量的方法，研究了关节空间角位置数值解法，给出了冗余自由度机构连续跟踪期望空间姿态时的关节空间角位置计算方法。在典型运动输入条件下，分别恒定速比关节位置算法和加权逆速度积分方法获得的关节位置和姿态误差进行仿真。为验证算法有效性，以卫星天线姿态稳定系统为研究对象，仿真结果表明，恒定速比关节位置算法在运动学逆解过程中能够保证天线姿态误差Frobenius(弗罗贝尼乌斯)范数小于10-4。仿真结果表明：本文提出的方法能够连续给出冗余自由度机构运动过程中的关节空间角位置的希望值，实现了逆运动学角位置层面的数值解和理想的控制精度。

冗余自由度机构； 卫星天线； 加权伪逆； 逆运动学； 运动规划； 恒定转速比； 关节角位置； 姿态稳定； 姿态误差； Frobenius范数

船舶在海上航行或无人飞行器执行远距离任务时，与陆地或指挥中心的通信都需要依靠地球同步通信卫星。上下行信号的传送，需要天线接收主瓣指向卫星，同时天线的极化方向与卫星发射信号的极化方向需要保持一致[1]。因此天线姿态稳定系统被广泛使用，实现卫星天线与地球同步卫星之间的相对稳定。接收天线主瓣方向与极化方向之间是垂直的，利用这两个方向构成了一个来描述天线姿态的坐标系。理论上，坐标系在空间的定点转动有3个自由度，具有3个正交自由度的运动机构就可以实现天线姿态的伺服。现实中，为了克服运动机构奇异点问题，增加系统灵活性，采用具有4个正交运动轴的框架结构，本文研究对象是一个冗余自由度的天线姿态稳定系统。

联接天线与运载器的稳定系统将实现了天线相对于运载器的姿态运动，同时也建立了天线姿态与稳定系统关节空间的约束。对与任务空间而言使用n维向量X，X∈Rn描述规划的任务；对于关节空间，使用m维向量Θ，Θ∈Cm描述关节变量。当m>n时系统即为冗余自由度系统，可以实现关节的避障、克服奇异点、完成第二任务等要求，例如实现机器人更好的灵活性[2-3]，躲避障碍[4-7]、防止关节运动超界[8-11]、避免奇异性[12-15]、减小关节运动速度及优化动力学性能等[16-17]。对于本文中的姿态随动系统，当系统自由度为m=4时，形成姿态冗余自由度系统。对于冗余自由度机器人而言，运动学逆解有无穷多组。从无穷多组解中选择合适解的过程，也就是冗余自由度问题的求解过程。成果包括：几何求解方法[18-19]；Pieper证明了对于连续3个轴交于一点的冗余自由度机构，一定具有解析解[20]；Liegeois及梯度投影法把冗余自由度问题的逆解分为了特解和齐次解，利用关节的自运动实现冗余自由度机器人的次任务[21-23]；Whitney最早采用了这种广义逆形式来求取一般机器人的逆运动学及轨迹规划问题，获得了关节角速度向量的模最小解，也成为基于速度层面求解冗余自由度问题的基础[24]；Wampler最早提出了加权最小二乘法实现主任务与自运动的同步完成[25]；Lee和Badler提出了可以选择合适的关节角度作为参数，表示其余的角度，进行解析解的求解方法[26-27]。

以上方法中，除几何方法和文献外，主要集中在对关节空间的速度求解，对于关节空间的位置求解研究较少，在对于机构末端轨迹有高精度要求的应用场合，并不适用[19]。

本文提出了一种建立在加权伪逆速度解基础上的恒定速比关节位置算法，使用该方法直接计算得到关节空间中的角位置值，计算误差不累计，能够保证天线姿态在运动过程中达到期望精度。

1 姿态稳定系统的描述

运载器、天线、稳定系统和卫星之间相互姿态关系与各坐标系的定义如图1所示。图1中，建立当地地理坐标系为东北天坐标系XeYeZe；运载器坐标系为与其固联的前右上坐标系XbYbZb；与天线固连的坐标系为XaYaZa，其中Za为天线接收主波瓣方向，Ya为天线极化方向；建立与卫星固联的坐标系XsYsZs，其中Zs为卫星发射信号主瓣的相反方向，Ys为卫星发射天线极化方向。

天线姿态稳定问题等效于将天线本体坐标系XaYaZa稳定在地理坐标系XeYeZe中，当运载器坐标系XbYbZb发生变化时，通过稳定系统的伺服运动，保证坐标系XaYaZa与卫星固联坐标系XsYsZs的一致。

图1 运载器、天线及卫星坐标系相互关系描述图Fig.1 Description of the relationship between carrier, antenna and satellite coordinate system

设每次一转动的姿态转移矩阵分别为T1、T2、T3和T4，而每一个姿态转移矩阵又仅与一个对应的关节角度有关，得到运载器坐标系与天线本体坐标系之间的姿态转移矩阵为

bTa=T1·T2·T3·T4=f(θ1,θ2,θ3,θ4)

(1)

式中：θ1～θ4为姿态稳定系统的4个运动轴角度值。尝试建立图1中各坐标系之间的转换关系，这里用姿态转移矩阵eTb表示坐标系XbYbZb在坐标系XeYeZe中的表达。同理利用姿态转移矩阵写出各坐标系之间转换关系为

bTs=bTe·eTs=bTa

(2)

式中：eTs为卫星在地理坐标系下的表示，对于给定的地球同步卫星目标，在地球上某一点观察卫星，可以认为eTs保持不变；eTb为运载器在地理坐标系下的表示；bTa为天线在运载器坐标系下的表示；运载器运动，eTb发生变化，系统完成的功能就是调整θ1～θ4，更新bTa使得(2)式保持成立[1]

bTa=bTe·eTs=N(t)

(3)

式中：将卫星在运载器坐标系下的连续运动用矩阵序列N(t)表示。则本文问题转化为如何在如式(4)的约束条件下，进行冗余自由度机构逆运动学求解。

f(θ1,θ2,θ3,θ4)=N(t)

(4)

矩阵N(t)的值则可通过姿态的运动和卫星参数求得：

这里假设卫星方位角度为Ae，俯仰角度为Ee，极化角约定为绕卫星姿态矩阵Zs轴经过σe的转动后，使得Ys轴与地面平行，此时的σe为极化角，并约定极化角正方向为逆时针方向，且确定极化角时以坐标系的Ys轴正方向为准。ψ为运载器航向角，φ为运载器的纵摇角，γ为运载器的横摇角。

2 加权意义下的速度最优解

2.1雅克比矩阵的建立

天线姿态稳定是一个姿态跟踪控制过程。将初始条件位置时刻定义为t0时刻，此时的天线姿态矩阵为N(t0)。随着运载器运动，在t1时刻运载器坐标系下的天线姿态应更新至N(t1)，在任意时刻姿态矩阵为N(tn)。根据机构的约束关系有：

f(Θn)=N(tn)

(5)

Θn=f-1(N(tn))

(6)

对于式(6)所描述的冗余自由度运动学逆解，f-1(·)代表函数f(·)的反函数。利用系统的Moore-Penrose的广义逆[24](简称M-P逆)可求得速度层面的最优解，本文中提出的恒定速比关节位置算法也是从伪逆速度的角度出发的。首先建立由雅克比矩阵J描述的关节运动与天线姿态运动的关系为

(7)

2.2天线坐标系的转动速度

在前节中已经获得了卫星姿态的描述矩阵eTs，且通过惯性测量单元得到了t时刻的运载器姿态的欧拉角，求得了eTb。利用式(2)直接得到连续变化的矩阵序列N(t)。对于任意相连时刻的N(tn)和N(tn+1)=N(tn+Δt)有：

N(tn+Δt)=Rot(q,dθ)N(tn)

(8)

式中：

Rot(q,dθ)=Rot(x,δx)Rot(y,δy)Rot(z,δz)

(9)

Rot(q,dθ)代表坐标系绕某一个轴q转动了dθ角度，Δt为两次测量测量之间的时间间隔；Rot(x,δx)为绕当前坐标系的x轴转动δx的微分运动，相应地Rot(y,δy)和Rot(z,δz)为绕另外两个轴的微分运动。当δx、δy和δz较小时，忽略二阶小量可得

(10)

2.3瞬时关节最优速度求解

许多文献针对冗余自由度机构在速度层面的解进行了深入的探讨，其中具有加权M-P逆解具有各关节同步运动，无关节自运动的特点，本文选用此方法作为关节速度的逆解方法[25]。对于四轴稳定系统，第1～4轴所消耗的能量逐渐减小，按照权值系数逐步降低的原则为转动关节分配转动权值，建立权值向量

(11)

此时的关节速度可表示为

(12)

3 恒定速比位置更新方法

3.1姿态跟踪系统的特点

姿态跟踪系统处在另外一种情况，一方面对于末端姿态精度有较高要求；另一方面，在给定一个期望姿态的情况下，希望系统够以最快的速度实现姿态更新，实现对目标位置。

对于姿态跟踪系统，姿态测量单元每秒更新20～50次姿态信息，每个更新的姿态信息对应姿态矩阵序列N(t)中的一个矩阵。每当新的姿态矩阵出现，要求跟踪系统能够尽可能快速地实现天线姿态与期望姿态N(t)的一致。对于两个姿态矩阵之间的姿态转移过程则不关注，不需要对两个姿态之间的轨迹做精确规划，而是将关注的重点放在对应每一个期望姿态所对应的关节角上，即对式(6)快速求解的方法，获得运动学逆解。

3.2现有方法的不足

2)姿态跟踪系统的控制不同于离线规划或示教方式规划的机器人系统，后者可以根据精度要求详细地规划轨迹上的每一个点。姿态跟踪系统关注的是期望姿态的快速建立和姿态建立精度，而对于给定的期望姿态如何直接求解关节的角位置尚没有成熟的方法。虽然可以建立姿态误差的反馈，多次计算伪逆，直至末端位姿满足精度要求。但计算量大，收敛速度慢，不适合进行在线规划。

3.3恒定速比位置算法

本文提出了一种基于关节最优速度比恒定和姿态约束方程的关节位置算法，简称为恒定速比(CRR)位置算法。

对于给定末端坐标系运动速度，基于式(12)可获得关节瞬时最优速度。冗余自由度最优速度意义本质上是将末端姿态运动合理地分解到各运动关节。角位移是有角速率和运动时间决定的，在相同的采样时间间隔上，如果角速率保持不变，那么关节空间角位移将正比于角速率。各轴的瞬时速率值，也代表了各轴瞬时角度运动量的大小。各关节瞬时速度之比，即瞬时关节运动量之比。将姿态转移过程中关节运动速度比设定为恒定，即恒定速比，是本文算法的核心思想。将每一个关节的运动速度通过某一个关节运动速度和关节之间速度的比值关系表示，并设定这一速比关系在每一次姿态转移过程中保持恒定。通过这样的表达，可以将采样周期内关节角位置增量表达出来，将这个待求角度增量的表达式直接带入姿态约束方程，求得每个关节角的增量值，进而获得关节角位置。

将此方法应用于天线姿态稳定系统，求解关节角位置的步骤为：

1)根据卫星的参数，计算矩阵eTs；

2)与运载器固连的姿态测量单元检测到运载器的运动姿态，利用欧拉角ψ、φ和γ描述航向角、纵摇角和横摇角，计算得到bTe，结合eTs求得此时卫星在运载器坐标系下的姿态矩阵记为N(tn+1)，相应的前一时刻的姿态矩阵记为N(tn);

3)在工作空间中进行运动微分的提取。计算连续两个位姿转动之间的差值，即两位姿之间的姿态转移矩阵为

N(tn+1)=Rot(q,dθ)N(tn)

(13)

(14)

5)根据建立的关节空间变量与工作空间之间的关系，在第n次姿态测量的时刻，有：

f(Θn)=N(tn)

(15)

那么对于tn+1有：

f(Θn+1)=f(Θn+ΔΘ)=

f(Θn+Φ·Δt)=N(tn+1)

(16)

6)建立恒等速比的关系式为

(17)

利用第一轴表示其余轴的速度，有：

(18)

式中:k2、k3、k4为关节速度比例系数，在单次关节位置计算过程中保持不变，即保持恒定速比。关节增量表达式可写为

(19)

7)将式(18)的增量表达式代入式(15)的方程中，将Θ作为初始，利用牛顿迭代方法即可获得ξ1。将ξ1代入式(19)，求得ΔΘ。

将算法总结为流程图，如图2所示。

图2 恒定速比位置更新算法流程图Fig.2 Flow chart of CRR position update algorithm

4 姿态稳定系统运动仿真结果

4.1仿真条件

本文为了验证算法在运载器大动态特性下的结果，设定了一个跟踪目标，并仿真运载器的运动规律，观察关节空间的运动情况及末端位姿的误差。仿真时间间隔为ΔT=0.01 s。

对跟踪目标做如下假设：目标方位角度：223.49°，俯仰角：60°，极化角度：-28.8°。按照本文中建模规则，初始位置下目标对应的关节角度为

(20)

运载器的航向角、纵横摇角运动组照如下的描述为

ψ=40sin[(2π/30 000)t]+0.01t

φ=40sin[(2π/800)t]

γ=45sin[(2π/1 000)t]

(21)

4.2伪逆方法仿真结果

按照仿真条件，分别使用恒等速比方法和加权M-P逆方法进行了稳定系统关节角、关节角速度和姿态误差的计算。

直接利用关节角速度积分或的关节角位置，如图4所示。

将计算得到的各关节角度带入运动学方程计算，获得末端的姿态矩阵，将此姿态矩阵与期望姿态矩阵进行做差，计算差矩阵的Frobenius范数值，如图5所示。

可以看到单纯利用伪逆速度积分获得的末端姿态矩阵误差会快速发散。

使用恒定速比方法进行计算，计算结果为每一步的关节角更新位置，直接获得关节角位置变化曲线如图6所示。

图3 根据M-P逆求得的各轴速度值Fig.3 The solution of joints velocity based on the M-P inverse

图4 对M-P逆生成速度进行积分得到的角度值Fig.4 The joints values obtained by using M-P inverse velocity integral

图5 姿态误差矩阵Frobenius范数值Fig.5 Frobenius norm values of attitude error matrix

本文解决的重点在关节空间角位置的计算，暂时不考虑关节的限位问题。在恒定速比方法求解过程中每一步计算得到的M-P逆速度如图7所示。

同样地根据各关节角度值进行运动学正计算，将获得姿态与期望姿态进行做差，求得误差矩阵的Frobenius范数曲线如图8所示。由于是利用恒速比方法直接求得的位置值，关节角的速度值是通过位置差分得到的平均速度值，差分结果如图9所示。

图6 利用本文的恒速比方法获得角度值Fig.6 The joints values obtained by using RCC approach

图7 恒速比方法中的各关节的M-P逆转速值Fig.7 The M-P inverse velocity of joints in RCC approach

图8 姿态误差矩阵Frobenius范数值Fig.8 Frobenius normal values of attitude error matrix

图9 恒定速比方法中位置差分得到的速度值Fig.9 The velocity of joints obtained by positiondifference in RCC approach

从图9误差结果可以看出，利用恒定速比方法计算关节角位置，机构的末端姿态误差不累加。虽然在某些位置出现误差增加点，但整体末端姿态始终保持较高精度。同时，姿态精度取决于算法中牛顿迭代的精度，换而言之，根据末端姿态精度要求的不同，可以灵活选择牛顿迭代法的精度。

5 结论

1)建立了具有4个运动自由度的天线稳定系统运动模型，给出了天线运动、运载器运动及目标卫星之间的关系；给出了天线姿态矩阵序列的计算方法。

2)分析了现有的冗余自由度运动学逆解方法在轨迹跟踪方面的不足。提出了一种基于关节恒定速比的关节位置算法，给出了算法的流程图。

3)对基于M-P逆的速度求解方法和本文提出的方法进行了仿真验证，仿真结果表明恒定速比方法得到的关节角，保证了系统末端姿态在运动过程中具有足够的精度。

[1] 吴鹏. 船用卫星天线姿态稳定系统关键技术研究[D].哈尔滨：哈尔滨工程大学, 2012： 17-23.

WU Peng. Key technologies of marine satellite antenna attitude stabilization system[D]. Harbin: Harbin Engineering University, 2012: 17-23.

[2] GUO Dongsheng, ZHANG Yunong. Acceleration-level inequality-based MAN scheme for obstacle avoidance of redundant robot manipulators[J]. IEEE transactions on industrial electronics, 2014, 61(12): 6903-6914.

[3] ZHANG Y, LI W, LIAO B,et al. Analysis and verification of repetitive motion planning and feedback control for omnidirectional mobile manipulator robotic systems[J]. Journal of intelligent & robotic systems, 2014,75 (3-4): 393-411.

[4] CAPISANI L M, FERRARA A. Trajectory planning and second-order sliding mode motion/interaction control for robot manipulators in unknown environments[J]. IEEE trans ind electron, 2012, 59(8): 3189-3198.

[5] XIAO L, ZHANG Y. A new performance index for the repetitive motion of mobile manipulators[J]. IEEE trans cybern, 2014,44(2)： 280-292.

[6] FLACCO F, DE LUCA A, KHATIB O. Control of redundant robots under hard joint constraints: saturation in the null space[J]. IEEE transactions on robotics, 2015, 31(3): 637-654.

[7] CHEN Dechao, ZHANG Yunong. A hybrid multi-objective scheme applied to redundant robot manipulators[J]. IEEE transactions on automation science and engineering, 2017，14(3)：1337-1350.

[8] FLACCO F, DE LUCA A. Optimal redundancy resolution with task scaling under hard bounds in the robot joint space[C]//Proc IEEE Int Conf Robot Seoul Korea, 2013: 3969-3975.

[9] ZHANG Dongsheng, XU Yundou , YAO Jiantao,et al. Kinematics, dynamics and stiffness analysis of a novel 3-DOF kinematically/actuation redundant planar parallel mechanism [J]. Mechanism and machine theory, 2017, 116:203-219.

[10] KUMAR A, KUMAR V. Evolving an interval type-2 fuzzy PID controller for the redundant robotic manipulator [J]. Expert systems with applications , 2017,73:161-177.

[11] KIM J, MARANI G, CHUNG W K, et al. A general singularity avoidance framework for robot manipulators: task reconstruction method[C]//IEEE International Conference on Robotics and Automation, New Orleas, 2004, 5: 4809-4814.

[12] MARANI G, KIM J, YUH J, et al. A real-time approach for singularity avoidance in resolved motion rate control of robotic manipulators[C]//IEEE International Conference on Robotics and Automation, Washington DC, 2002, 2: 1973-1978.

[13] ZORJAN M, HUGEL V. Generalized humanoid leg inverse kinematics to deal with singularities[C]∥IEEE International Conference on Robotics and Automation Karlsruhe Germany, 2013: 4791-4796.

[14] JAKUBIAK J, RATAJCZAK A. Kinematics and motion planning of the multi-bar system[J]. Journal of intelligent and robotic systems, 2014, 75(1): 5-15.

[15] RANDY C H, RODNEY G R, ANTHONY A M, et al. Designing a failure-tolerant workspace for kinematically redundant robots[J]. IEEE transactions on automation science and engineering, 2015, 12(4): 1421-1432.

[16] ZHANG Yunong, ZHANG Yinyan, CHEN Dechao, et al. From davidenko method to zhang dynamics for nonlinear equation systems solving[C]∥IEEE Transactions on Systems, Man, and Cybernetics: Systems. Budapest Hungary, 2016.

[17] LI S, LI Y. Nonlinearly activated neural network for solving timevarying complex Sylvester equation[J]. IEEE trans. cybern, 2014, 44(8): 1397-1407.

[18] ASFOUR T, DILLMANN R. Human-like motion of a humanoid robot arm based on a closed-form solution of the inverse kinematics problem[C]∥IEEE/RSJ International Conference on Intelligent Robots and Systems. Las Vegas Nevada, 2003, 2: 1407-1412.

[19] 王英石.冗余机器人的运动学及轨迹规划的研究[D].天津：南开大学, 2014.

WANG Yingshi. Research on the kinematics and trajectory planning of redundant robots[D].Tianjin: Nankai University, 2014.

[20] PEIPER D L. The kinematics of manipulators under computer control[R]. Stanford Univ Calif Dept of Computer Science, 1968: 32-66.

[21] LIEGEOIS A. Automatic supervisory control of the configuration and behavior of mutibody mechanisms[J]. IEEE trans sys man cyber,1977, SMC-7(12): 868-871.

[22] JIN Hongzhe, ZHANG Hui, LIU Zhangxing,et al. A synthetic algorithm for tracking a moving object in a multiple-dynamic obstacles environment based on kinematically planar redundant manipulators [J]. Mathematical problems in engineering, 2017(2): 1-15.

[23] DUBEY R V, EULER J A, BABCOCK S M. An efficient gradient projection optimization scheme for a seven-degree-of-freedom redundant robot with spherical wrist[C]∥IEEE international conference on robotics and automation, philadelphia, PA, USA, 1988: 28-36.

[24] WHITNEY D E. The mathematics of coordinated control of prosthetic arms and manipulators[J]. ASME journal of dynamic systems, measurement, and control, 1972, 94(4): 303-309.

[25] WAMPLER C W. Manipulator inverse kinematic solutions based on formulations and damped least-squares methods[J]. IEEE transactions on systems man and cybernetics, 1986, 16(1): 93-101.

[26] CHEN Yujen, JU Mingyi, HWANG Kaoshing. A Virtual torque-based approach to kinematic control of redundant manipulators[J]. IEEE Transactions on industrial electronics, 2017, 64(2): 1728-1736.

[27] BADLER N I, TOLANI D. Real-time inverse kinematics of the human arm[J]. Center for human modeling and simulation, 1996, 73: 392-401.

本文引用格式：吴鹏，洪娟，陈广，等. 串联冗余自由度机构关节角位置的规划算法[J]. 哈尔滨工程大学学报， 2017， 38(10): 1623 -1629，1641.

WU Peng, HONG Juan, CHEN Guang, et al. Joints position planning approach for series redundant DOF mechanisms[J]. Journal of Harbin Engineering University， 2017， 38(10): 1623 -1629，1641.

JointanglepositionplanningapproachforseriesredundantDOFmechanisms

WU Peng1, HONG Juan2, CHEN Guang1, XU Bo1

(1.College of Automation, Harbin Engineering University, Harbin 150001, China; 2.Beijing Institute of Aerospace Control Device, Beijing 100854, China)

This paper proposes a method that can solve the joint angle position of a mechanism with series redundant DOFs in real time, solve the problem of attaining angular position in the inverse kinematics of the mechanism with series redundancy, and provide the numerical solution of inverse kinematics. According to weighted pseudo-inverse theory, this paper adopted the instantaneous optimum speed to represent the angle increment of joint space, studied the numerical solution of the joint space position, and presented the computing method for the joint space's angle position as the redundant DOF mechanism tracks the desired space attitude continuously. Under typical motion input conditions, the joint position and attitude errors obtained by constant speed ratio joint position algorithm and weighted inverse velocity integral method were simulated. To verify the effectiveness of the algorithm, the satellite antenna attitude stabilization system was taken as the research target. Simulation results show that the constant speed ratio joint position algorithm can ensure that the attitude error of the antenna Fresenius norm is less than 10-4in the process of inverse kinematic solution. Simulation results show that the proposed method can give the desired value of joint space position continuously in the movement process of a mechanism with redundant DOFs. The numerical solution of the position of the inverse kinematics can be realized, and control precision is ideal.

redundant manipulators; satellite antenna; weighted pseudo inverse; inverse kinematics; motion planning; constant rotation ratio (CRR); joint position; attitude stabilization; attitude error; Frobenius norm

10.11990/jheu.201607014

http://www.cnki.net/kcms/detail/23.1390.U.20170816.1550.042.html

TP241.3

A

1006-7043(2017)10-1623-08

2016-07-05. < class="emphasis_bold">网络出版日期

日期：2016-08-16.

国家自然科学基金项目(61203225，61633008)；黑龙江省自然科学基金项目(QC2014C0669)；中央高校专项基金项目(HEUCF041702)；海洋工程国家重点实验室开放课题项目(1616).

吴鹏(1978-)，男，讲师，博士；徐博(1984-)，男，副教授，研究生.

徐博，E-mail:xubocarter@sina.com.