基于线调频小波路径追踪算法与EEMD的齿轮箱复合故障诊断方法

李 蓉，于德介，陈向民，刘 坚

（湖南大学汽车车身先进设计制造国家重点实验室，长沙 410082）

基于线调频小波路径追踪算法与EEMD的齿轮箱复合故障诊断方法

李 蓉，于德介，陈向民，刘 坚

（湖南大学汽车车身先进设计制造国家重点实验室，长沙 410082）

针对变转速下的齿轮箱中复合故障的故障特征提取，提出一种基于线调频小波路径追踪算法与集合经验模式分解的齿轮箱复合故障诊断方法。该方法先用线调频小波路径追踪算法从原始振动信号中提取转频曲线，根据转频曲线对原始振动信号进行等角度重采样，将时域信号转化为角域信号，再对角域重采样信号进行集合经验模式分解，根据相关系数选取合适的内禀模态函数，最后对所选取的内禀模态函数分量进行Hilbert包络谱分析，根据包络谱进行齿轮箱复合故障诊断。通过算法仿真和应用实例对包含齿轮局部故障和轴承局部故障的变转速齿轮箱复合故障进行分析，结果表明，该方法在无转速计的情况下能有效地提取变转速齿轮箱复合故障的特征。

线调频小波；阶次跟踪；集合经验模式分解；齿轮箱；复合故障

齿轮箱是机械系统中的重要部件，因其工作环境恶劣，且常处于变载荷下运行，容易因疲劳磨损而发生局部故障。在实际生产中，故障往往不是单独出现的，某些故障常常会诱发其他故障的发生［1］，因此，对齿轮箱的复合故障进行研究具有重要的实际意义。

齿轮、轴承作为齿轮箱中的主要部件，在齿轮箱的机械传动中起着重要作用。当齿轮出现断齿、裂纹等局部故障时，信号中会出现转频调制现象［2］；而当轴承内圈、外圈出现局部故障时，会出现周期性的冲击，冲击出现的频率（调制频率）为轴承内圈、外圈的通过频率［3］。当齿轮箱中同时出现齿轮局部故障和轴承局部故障时，由于齿轮故障特征调制频率与轴承故障特征调制频率不同，因此，根据调制频率的差异，即可实现对包含齿轮局部故障和轴承局部故障的齿轮箱复合故障进行诊断。

Hilbert解调为常用的包络分析方法，能有效的提取信号中的调制信息。但由于齿轮箱的工作环境往往比较恶劣，其故障特征信息往往淹没在在强背景噪声环境中，不适合直接进行Hilbert分析，需预先进行降噪处理。集合经验模态分解（Ensemble Empirical ModeDecomposition，EEMD）利用了高斯白噪声具有频率均匀分布的统计特性，在信号中加入了高斯白噪声，使得信号在不同尺度上具有连续性，避免了经验模态分解（Empirical Mode Decomposition，EMD）方法中由于内禀模态函数（Intrinsic Mode Function，IMF）的不连续性而造成的模式混淆现象［4］，且在抗噪能力上，EEMD要优于小波和EMD［5－6］。但由于EEMD继承了EMD的分解特性，即按照频率从高到低进行分解，当信号为频率大范围波动的非平稳信号时，直接进行EEMD分解也容易造成模式混淆，因而，需预先对信号进行平稳化处理。

工程实际中常采用阶比跟踪方法［7－9］来实现信号的平稳化，其基本思想是将时域非平稳振动信号进行等角度重采样转化为角域平稳振动信号。常用的阶比跟踪方法有硬件阶比跟踪法、计算阶比跟踪法［10－11］、基于瞬时频率估计的阶比跟踪法［12］等。其中，基于瞬时频率估计的阶比跟踪法因无需安装角度编码盘、转速计等硬件近年来得到了广泛的关注。而在基于瞬时频率估计的阶比跟踪法中，其核心问题是如何从原始振动信号中准确地提取出转速信号。Emmanue等［13］近年提出了线调频小波路径追踪算法，该算法通过对线调频小波图中的线调频小波原子进行连接，自适应的获得频率呈曲线变化的信号分量。由于该方法具有精度高、抗噪能力强等优点，近年被引入机械故障诊断［14－15］，能自适应地从振动信号中提取出转频信息。

针对变转速下齿轮箱复合故障诊断问题，本文结合线调频小波路径追踪算法、阶次跟踪和EEMD方法，提出了基于线调频小波路径追踪与EEMD的齿轮箱复合故障诊断方法。该方法先采用线调频小波路径追踪算法自适应地从变转速齿轮箱原始振动信号中提取转速信号；然后利用该转速信号对原始振动信号进行等角度重采样，将时域信号转化为角域信号；再对角域重采样信号进行EEMD分解，根据相关系数选取合适的IMF分量，最后对所选取的IMF分量进行Hilbert包络谱分析，根据包络谱进行齿轮箱复合故障诊断。算法仿真和应用实例表明，本文方法在无转速计的情况下，能有效地提取出包含齿轮局部故障和轴承局部故障的齿轮箱复合故障的故障特征。

1 线调频小波路径追踪算法

线调频小波路径追踪算法采用的多尺度线调频基元函数库如下［13］：

式中：D为基元函数库；

haμ，bμ，I（t）为多尺度线调频基元函数；

k为动态时间段序号，k＝0，1，…，2j－1；

N为分析信号的采样长度；

j为分析尺度系数，j＝0，1，…，log2（N－1）；

Kaμ，bμ，I为归一化系数，使得

aμ为频率偏置系数；

bμ为调频率；根据采样定理aμ＋2bμt应该小于fs／2；

1I（t）为矩形窗函数，当t∈I时为1，当tI时为0。

式（1）定义的多尺度线性调频基函数在动态分析时间段内的瞬时频率为aμ＋2bμt。通过多尺度线性调频基函数对信号进行逐段投影分析，计算获得每个时间分析段I内的最大投影系数和对应的线调频基元函数，该基元函数即为在时间分析段I中与分析信号最为相似的频率成分。

当信号与多尺度线性调频基函数越相似时，其投影系数也越大，基元函数的能量也越大，因此要求找到一种动态分析时间段连接方法，在该连接方法下整个分析时间内连接的所有基元函数信号的总能量最大，即：

（1）初始化。以i为时间支持区序号，di为第i个时间支持区之前分解信号的总能量，pi为连接到第i个时间支持区的前置时间支持区序号，ei为第i个时间支持区最大投影系数对应的分解信号的能量，初始化时，置di＝0，pi＝0；

（2）对于动态分析时间段集合｛Ii，i∈Z｝中的每一个元素Ii，查找出与其相邻的所有下一个动态分析时间段集合｛Ij｝，即｛Ij｝中所有元素的起始时间与Ii相邻。如果

di＋ei＞dj（4）

有：

dj＝dj＋eipj＝i（5）

2 EEMD原理

EMD方法是一种自适应信号分解方法，适合于处理非线性、非平稳信号。它通过一种“筛”的过程从被分析信号中提取其本身固有的一族本征模态函数（IMF）［16］。

EEMD方法是EMD方法的改进，它利用了高斯白噪声具有频率均匀分布的统计特性，向信号中加入高斯白噪声，使得信号在不同尺度上具有连续性，这样可以避免EMD方法中由于IMF的不连续性而造成的模式混淆现象。EEMD的分解过程如下：

（1）在原始信号中加入高斯白噪声；

（2）将加噪后的信号进行EMD分解，得到一族IMF；

（3）每次加入相同幅值的不同高斯白噪声，重复步骤（1）和步骤（2）；

（4）把分解得到的各个IMF的均值作为最终的结果，即：式中：M为加入高斯白噪声的次数，cij（t）为第i次加噪后经EMD分解得到的第j个IMF分量，cj（t）为对原始信号进行EEMD分解得到的第j个分量。

3 基于线调频小波路径追踪与EEMD的齿轮箱复合故障诊断

由于线调频小波路径追踪算法能自适应地从原始振动信号中精确地估计出信号中的转速信息，不需要安装任何与转速相关的硬件；而EEMD方法解决了EMD中的模态混叠现象，能将故障成分从信号中有效分离出来，增强信噪比，因此，本文将线调频小波路径追踪、阶次跟踪和EEMD方法相结合，提出了基于线调频小波路径追踪与EEMD的齿轮箱复合故障诊断方法。该方法先采用线调频小波路径追踪算法从变转速齿轮箱原始振动信号中估计转速信号；然后根据该转速信号对原始振动信号进行等角度重采样，将时域信号转化为角域信号；再对角域重采样信号进行EEMD分解，同时，根据相关系数选取合适的IMF分量；最后对所选取的IMF分量进行Hilbert包络谱分析，并根据包络谱进行齿轮箱复合故障诊断。本文方法主要步骤如下：

（1）利用线调频小波路径追踪算法从原始振动信号x（t）中估计出转速曲线fr；

（2）根据估计出的转速曲线fr对原始振动信号x（t）进行等角度重采样，得到角域重采样信号；

（3）确定加入高斯白噪声的幅值k和总体平均次数M。一般情况下，k取0.2［4］，而M的取值理论上越大越好，但随着M的增大，计算量随之增加，本文中k取0.2，M取200。对角域重采样信号进行EEMD分解，得到各IMF分量；

（4）利用式（7）求取各IMF分量与角域重采样信号的相关系数；

（5）选择相关系数阈值σ（本文取0.5），对相关系数大于σ的IMF进行Hilbert包络谱分析，根据包络谱诊断齿轮箱复合故障。

4 算法仿真

为验证本文方法对变转速状况下齿轮箱复合故障诊断的有效性，构造包含齿轮局部故障和轴承局部故障的变转速齿轮箱复合故障仿真信号，如式（8）所示。式（8）中，x1（t）为被1倍转频调制的调幅调频信号，以模拟变转速下的齿轮故障信号，模拟齿数为25，其信号波形图如图1（a）所示；x2（t）为周期性的瞬态冲击信号，以模拟变转速下的轴承故障模拟信号，y（t）为单瞬态冲击信号，瞬态冲击的中心频率为1 200 Hz，衰减系数为－800；u（t）为单位阶跃函数；P表示xc（t）的最大值（即峰值）集合，Max（·）表示取最大值；tpi表示第i个峰值所对应的时间点；xc（t）为轴承故障调制信号，轴承故障阶次为f0＝3.4，即轴承故障特征频率为转频fr的3.4倍；fr为轴的转频，则x2（t）的信号波形图如图1（b）所示。n（t）为－4 dB的高斯白噪声，由MATLAB中的函数awgn产生，以模拟随机干扰。信号采样频率4 096 Hz，采样时长为1 s。式（8）合成信号的时域波形如图1（c）所示。图2为图1（c）合成信号的幅值谱，图中出现了频率“模糊”现象，故无法判断故障类型。

fr＝20＋2×π×sin（1.4×π×t）（14）

利用线调频小波路径追踪算法从图1（c）所示合成信号中估计瞬时转频，获取的估计瞬时转频曲线如图3中实线所示，而虚线则为实际瞬时转频。从图3中可看出，实线与虚线基本重合，说明线调频小波路径追踪算法能很好地从信号中提取转频曲线，具有较高的精度。

利用图3中的估计瞬时转频曲线对图1（c）所示信号进行角域重采样，得到的角域重采样信号如图4所示。

对角域重采样信号进行EEMD分解，得到各IMF分量，求取各IMF分量与角域重采样信号的相关系数，其中，IMF1～IMF5分量的相关系数如表1所示。

表1 角域重采样信号与各IMF的相关系数Tab.1 The cor relation coefficient of the angle domain resampling signal and each IMF

由于IMF1和IMF2与角域重采样信号较为相关（即大于阈值系数），故取IMF1与IMF2进行Hilbert分析，得到的包络谱图5所示。图5（a）中，在轴承故障特征阶次O0和二倍特征阶次2O0处存在明显峰值；图5（b）中，在转频阶次Or处峰值明显，即信号中出现了1倍转频调制现象，与实验设置相符，验证了本文方法的有效性。

图1 模拟齿轮故障信号、模拟轴承故障信号与合成信号Fig.1 The simulated signal of the fault gear、the simulated signal of the fault bearing and the composite signal

图2 合成信号的幅值谱Fig.2 The amplitude spectrum of the composite signal

图3 实际瞬态转频与估计瞬态转频对比Fig.3 The comparison between the actual rotate speed and the estimated rotate speed

图4 角域重采样信号Fig.4 The angle domain resampling signal

图5 IMF分量的包络谱Fig.5 The envelope spectrum of IMFs

图6 角域循环平稳信号的包络谱Fig.6 The envelope spectrum of the angle domain cyclostationary signal

对图4所示角域重采样信号直接进行Hilbert包络谱分析，得到的包络谱图如图6所示。图6中，在转频阶次Or处峰值突出，但在轴承故障阶次处无明显峰值，即无法从图6诊断出轴承故障。

5 应用实例

当齿轮箱中的轴承内圈或外圈出现裂纹、点蚀等局部故障时，会产生周期性的冲击信号，其冲击出现的频率（调制频率）为内圈或外圈的通过频率，内圈、外圈的通过频率计算公式分别如式（15）、式（16）所示［2］：式中：Z为滚动体个数，fr为轴的转频，d0为滚动体直径，D为节径，α为接触角。而当齿轮出现断齿、裂纹等局部故障时，信号中会产生转频调制现象。为验证本文方法提取齿轮箱复合故障的故障特征的有效性，在单级传动齿轮箱试验台上进行试验，试验台其简图如图7。试验齿轮为正齿轮，主动轴与从动轴齿数均为37。试验轴承1～4均为SKF 6307－2RS深沟球轴承，其内径d＝35 mm，外径D＝80 mm，滚动体个数Z＝8，接触角α＝0，滚动体直径d0＝13.5 mm，经计算，内圈、外圈对应的故障特征阶次分别为4.94、3.06，即内圈、外圈的故障特征频率分别为轴转频的4.94倍和3.06倍。为模拟齿轮箱齿轮、轴承复合故障，设置表2中齿轮箱复合故障类型进行试验。

表2 齿轮箱复合故障类型Tab.2 Compound fault type of gearbox

在齿轮2上整体切割掉一个齿，以模拟齿轮断齿故障，并用激光在齿轮2齿根处切割宽0.15 mm、深1 mm的槽，以模拟齿轮裂纹故障；同时用激光在轴承4的内圈和外圈上切割宽0.15 mm，深0.13 mm的槽，以模拟轴承内圈和外圈故障。为减少传递路径的影响，振动加速度传感器置于轴承4的轴承盖上，测取径向垂直方向上的振动加速度信号。试验用LMS数据采集设备采集振动加速度信号，采样频率为8 192 Hz，采样时长为4 s。

图7 试验台简图Fig.7 The diagram of the test rig

5.1 复合故障1振动信号分析

图8为存在齿轮2断齿故障和轴承4外圈故障的时域波形图，图中可看出，信号中存在明显冲击，但冲击之间的时间间隔逐渐变大，且幅值随时间变小，说明齿轮箱处于降速阶段。

利用本文方法对图8所示信号进行分析，得到的前5个IMF分量与角域重采样信号的相关系数分别为：0.773、0.446、0.276、0.194、0.145，故取IMF1进行Hilbert分析，得到的包络谱如图9所示。图9中，在转频阶次Or和二倍转频阶次2Or处出现峰值，说明齿轮箱出现转频调制现象，与齿轮断齿故障相符；同时在轴承外圈故障特征阶次及倍频阶次O0～5O0处出现明显峰值，表示轴承外圈出现局部故障。

5.2 复合故障2振动信号分析

图10为存在齿轮2断齿故障和轴承4内圈故障的时域波形图，图示信号中存在冲击。

利用本文方法对图10所示信号进行分析，得到的前5个IMF分量与角域重采样信号的相关系数分别为：0.759、0.345、0.232、0.170、0.122，故取IMF1进行Hilbert分析，得到的包络谱如图11所示。图11中，在转频阶次及倍频阶次Or～3Or处出现峰值，说明齿轮箱出现转频调制现象，与齿轮断齿故障相符；同时在轴承外圈故障特征阶次Oi和二倍外圈特征阶次2Oi处出现明显峰值，表示轴承内圈出现局部故障。

图8 齿轮断齿与轴承外圈故障振动信号Fig.8 The vibration signal of the test rig with broken gear and outer race cracked bearing

图9 IMF1的包络谱（齿轮断齿与轴承外圈故障）Fig.9 The envelope spectrum of IMF1（broken gear and outer race cracked rolling bearing）

图10 齿轮断齿与轴承内圈故障振动信号Fig.10 The vibration signal of the test rig with broken gear and inner race cracked rolling bearing

5.3 复合故障3振动信号分析

图12为存在齿轮2裂纹故障和轴承4外圈故障的时域波形图，图中信号中存在明显冲击，冲击的幅值由大变小，冲击之间的时间间隔由小变大，说明信号处于降速阶段。

图11 IMF1的包络谱（齿轮断齿与轴承内圈故障）Fig.11 The envelope spectrum of IMF1（broken gear and inner race cracked rolling bearing）

图12 齿轮裂纹与轴承外圈故障振动信号Fig.12 The vibration signal of the test rigwith cracked gear and outer race cracked rolling bearing

图13 IMF1的包络谱（齿轮裂纹与轴承外圈故障）Fig.13 The envelope spectrum of IMF1（cracked gear and outer race cracked rolling bearing）

利用本文方法对图12所示信号进行分析，得到的前5个IMF分量与角域重采样信号的相关系数分别为：0.746、0.485、0.289、0.211、0.131，故取IMF1进行Hilbert分析，得到的包络谱如图13所示。图13中，在转频阶次Or、二倍转频阶次2Or和四倍转频阶次4Or处出现峰值，说明齿轮箱出现转频调制现象，与齿轮裂纹故障相符；同时，在轴承外圈故障特征阶次及倍频阶次O0～6O0处出现明显峰值，表示轴承外圈出现局部故障。

5.4 复合故障4振动信号分析

图14为存在齿轮2裂纹故障和轴承4内圈故障的时域波形图，图示信号中存在冲击。

图14 齿轮裂纹与轴承外圈故障振动信号Fig.14 The vibration signal of the test rig with cracked gear and inner race cracked rolling bearing

图15 IMF1的包络谱（齿轮裂纹与轴承内圈故障）Fig.15 The envelope spectrum of IMF1（cracked gear and inner race cracked rolling bearing）

利用本文方法对图14所示信号进行分析，得到的前5个IMF分量与角域重采样信号的相关系数分别为：0.759、0.477、0.304、0.192、0.133，故取IMF1进行Hilbert分析，得到的包络谱如图15所示。图15中，在转频阶次及倍频阶次Or～4Or，6Or处出现峰值，说明齿轮箱出现转频调制现象，与齿轮裂纹故障相符；同时，在轴承外圈故障特征阶次Oi处出现明显峰值，表示轴承内圈出现局部故障。

6 结 论

（1）齿轮出现局部故障时振动信号中往往出现转频调制现象，而轴承出现局部故障时其振动信号则会被轴承通过频率调制，因而，可根据调制频率的不同，即可实现对包含齿轮局部故障和轴承局部故障的齿轮箱复合故障进行诊断。

（2）通过算法仿真和应用实例对包含齿轮局部故障与轴承局部故障的变转速齿轮箱复合故障进行分析，结果表明，本文方法可在无转速计的情况下，有效地提取齿轮箱中复合故障的故障特征。

［1］张含蕾，周洁敏，李 刚.基于小波分析的感应电动机复合故障诊断［J］.中国电机工程学报，2006，26（8）：159－162.

ZHANG Han-lei，ZHOU Jie-min，LI Gang.Mixed fault diagnosis based on wavelet analysis in induction motors［J］.Proceedings of the CSEE，2006，26（8）：159－162.

［2］丁 康，李巍华，朱小勇.齿轮及齿轮箱故障诊断实用技术［M］.北京：机械工业出版社，2005.

［3］Ho D，Randall R B R.Optimization of bearing diagnostic techniques using simulated and actual bearing fault signals［J］.Mechanical Systems and Signal Processing，2000，14（5）：763－768.

［4］Wu Z H，Huang N E，Ensemble empirical mode decomposition：a noise assisted data analysis method［J］.Advances in Adaptive Data Analysis，2009，1（1）：1－41.

［5］曹冲锋，杨世锡，杨将新.大型旋转机械非平稳振动信号的EEMD降噪方法［J］.振动与冲击，2009，28（9）：33－38.

CAO Chong-feng，YANG Shi-xi，YANG Jiang-xin.The EEMD de-noisingmethod in non-stationary vibration signal of large rotating machinery［J］.Journal of Vibration and Shock，2009，28（9）：33－38.

［6］陈仁祥，汤宝平，马婧华.基于EEMD的振动信号自适应降噪方法［J］.振动与冲击，2012，31（15）：82－86.

CHEN Ren-xiang，TANG Bao-ping，MA Jing-hua.Adaptive de-noising method based on ensemble empirical mode decomposition for vibration signal［J］.Journal of Vibration and Shock，2012，31（15）：82－86.

［7］Borghesani P，Pennacchi P，Randall R B，etal.Order tracking for discrete-random separation in variable speed conditions［J］.Mechanical Systems and Signal Processing，2012，30：1－22.

［8］Shao H，Jin W，Qian S.Order tracking by discrete gabor expansion［J］.IEEE Transactions on Instrumentation and Measurement，2003，52（3）：754－761.

［9］赵晓平，侯荣涛.基于Viterbi算法的Gabor阶比跟踪技术［J］.机械工程学报，2009，45（11）：247－252.

ZHAO Xiao-ping，HOU Rong-tao.Gabor order tracking based on Viterbialgorithm［J］.Journal of Mechanical Engineering，2009，45（11）：247－252.

［10］Bossley K M，Mckendrick R J，Harris C J，et al.Hybrid computed order tracking［J］.Mechanical Systems and Signal Processing，1999，13（4）：627－641.

［11］Saavedra P N，Rodriguez C G.Accurate assessment of computed order tracking［J］.Shock and Vibration，2006，13（1）：13－32.

［12］郭 瑜，秦树人，汤宝平，等.基于瞬时频率估计的旋转机械阶比跟踪［J］.机械工程学报，2003，39（3）：32－36. GUO Yu，QIN Shu-ren，TANG Bao-ping，et al.Order tracking of rotating machinery based on instantaneous frequency estimation［J］.Chinese Journal of Mechanical Engineering，2003，39（3）：32－36.

［13］Candès E J，Charlton P R，Helgason H.Detecting highly oscillatory signals by chirplet path pursuit［J］.Applied and Computational Harmonic Analysis，2008，24（1）：14－40.

［14］Peng F Q，Yu D J，Luo J S.Sparse signal decomposition method based onmulti-scale chirplet and its application to the fault diagnosis of gearboxes［J］.Mechanical System and Signal Processing，2011，25（2）：549－557.

［15］陈向民，于德介，罗洁思.基于线调频小波路径追踪阶比循环平稳解调的齿轮故障诊断［J］.机械工程学报，2012，48（3）：95－101.

CHEN Xiang-min，YU De-jie，LUO Jie-si.Gear tooth fault diagnosis by using order cyclostationary demodulating approach based on chirplet path pursuit［J］.Journal of Mechanical Engineering，2012，48（3）：95－101.

［16］程军圣，于德介，杨 宇.基于EMD的能量算子解调方法及其在机械故障诊断中的应用［J］.机械工程学报，2004，40（8）：115－118.

CHENG Jun-sheng，YU De-jie，YANG Yu.Energy operator demodulation approach based on EMD and its application in mechanical fault diagnosis［J］.Chinese Journal of Mechanical Engineering，2004，40（8）：115－118.

A compound fault diagnosismethod for gearboxs based on chirplet path pursuit and EEMD

LIRong，YU De-jie，CHEN Xiang-min，LIU Jian
（State Key Laboratory of Advanced Design and Manufacturing for Vehicle Body，Hunan University，Changsha 410082，China）

Aiming at the problem of extracting fault characteristics from a gearbox with compound faults under a condition of changing rotating speed，a compound fault diagnosismethod for gearboxs based on chirplet path pursuit and ensemble empiricalmode decomposition（EEMD）was proposed.Using the proposed method，the rotating speed was estimated from an original vibration signalwith chirplet path pursuitalgorithm.In order to transform a time domain signal to an angular domain one，the even angle resampling was performed to the original signal according to the obtained rotating speed.The angular domain resampled signal was decomposed with EEMD，and the suitable intrinsic mode function（IMF）was selected according to the correlation coefficient.The selected IMF was analyzed with Hilbert envelope spectrum and the compound fault diagnosis of the gearbox was executed according to the envelope spectrum.The results obtained with algorithm simulations and application examples of a gearbox with compound faults including gear local faults and bearing local faults under the condition of changing rotating speed showed that the proposed method can effectively extract fault characteristics from a gearbox with compound faults.

chirplet；order tracking；EEMD；gearbox；compound fault

TH113.1；TH165.3

A

湖南省科技计划资助（2012SK3184）；广东省省部产学研结合项目（2009B090300312）；高等学校博士学科点专项科研基金（20090161110006）资助项目

2012－11－07 修改稿收到日期：2013－03－06

李 蓉女，硕士，副教授，1976年生

于德介男，教授，博士生导师，1957年生