CN105573290B - 卷烟工厂超高速条盒包装机多工况过程在线监测与故障诊断方法 - Google Patents

Abstract

本发明涉及卷烟包装过程中意大利G.D公司超高速条盒包装机BV的统计建模、在线监测与故障诊断技术。卷烟工厂超高速条盒包装机多工况过程在线监测与故障诊断方法，本发明基于滑动时间窗口离线计算稳定度因子识别稳定工况和过渡工况，采用自适应k‑means聚类方法对稳定工况进行划分形成若干个稳定工况数据簇，再利用PCA方法对每个稳定工况数据簇独立建立统计监测模型。在线监测时，根据当前滑动时间窗口内数据的稳定度因子判断工况类型，在稳定工况下利用距离最小的聚类所对应的PCA监测模型进行实时监测，任一统计量超限时利用贡献图方法确定引起故障的主要过程变量，最终为超高速条盒包装机BV的在线监测与故障诊断提供了一种有效可行的方法。

Description

卷烟工厂超高速条盒包装机多工况过程在线监测与故障诊断 方法

技术领域

本发明涉及卷烟包装过程中意大利G.D公司超高速条盒包装机BV的统计建模、在线监测与故障诊断技术。

背景技术

随着近年来烟草行业全球化、新的法律法规以及外部环境的重大变化，中国烟草面临愈来愈严峻的竞争压力和社会压力。在保证产品质量的前提下，有效提升设备的智能化水平和高效运行能力成为卷烟工厂关注的焦点。智能化，即生产设备具有感知、分析、推理和决策等功能；高效运行，即在规定时间内消耗更少的资源生产更多、更好的产品。卷烟工厂设备的智能化水平和高效运行能力与其自动化程度密不可分，设备状态的在线监测与故障诊断技术是提升智能化水平、保证高效运行的重要手段。

意大利G.D公司超高速包装设备是中国烟草总公司“十二五”期间重点引进的目前世界上最高速、最先进的卷烟包装设备，其设计生产能力为1000包/分钟香烟。超高速包装设备由小盒包装机H1000、小盒薄膜封装机W1000和条盒包装机BV三部分组成。条盒包装机BV是一台高速自动外条包装机，可以包装各种不同特性和硬度的条包，其运行过程由二五平包成型、条盒纸松卷、条盒纸上胶、条盒纸输送、条盒纸包裹、条盒纸折叠、透明纸松卷、透明纸分切、透明纸输送、透明纸包裹、透明纸折叠成型、透明纸热封、条烟输出等部分组成。条盒包装机BV是超高速包装设备的重要组成部分，也是生产现场故障较为集中的部位。由于G.D公司超高速包装设备底层PLC开放的程序源代码和库文件非常有限，生产现场主要通过每班一组的数据报表、轮保时的人工检测、设备故障的计划维修和事后抢修等方式保证条盒包装机BV的完好性。2013、2014年杭州卷烟厂G.D公司超高速包装设备的运行效率分别为70％和70.08％，其中条盒包装机BV的条盒纸和透明纸堵塞、透明纸折角等故障是影响停设备运行效率的主要因素之一。目前采用的条盒包装机BV状态监测方式已经不能满足卷烟生产的实际需要，亟需一种切实有效的条盒包装机BV在线监测与故障诊断方法。

为解决G.D公司超高速包装设备底层PLC参数数据难以获取的问题，杭州卷烟厂(朱立明,钱杰,王琪,李钰靓,李勇.一种超高速包装机数据采集系统:中国,CN204056360U[P].2014-12-31.)提出了一种基于Beckoff软PLC的超高速包装设备底层PLC数据采集方法，通过对设备底层PLC源代码和库文件的分析，能够采集条盒包装机BV的部分参数数据，从而为基于多元统计分析方法实施条盒包装机BV统计建模、在线监测与故障诊断提供了丰富的运行数据。

条盒包装机BV生产过程中，操作人员根据目测的设备运行状态和产品外观质量频繁调整BV车速设定值，使得生产过程具有多个运行工况。由于条盒包装机BV的高速运行，从一个稳定运行工况到另外一个稳定运行工况的时间非常短(如BV车速从89条/分钟降低到0条/分钟仅需10秒，即在一个采样周期内完成)，因此本发明仅考虑各个不同的稳定运行工况。不同稳定运行工况之间数据的均值、方差、相关关系等特征变量会有明显的变化，直接应用PCA方法进行统计建模和在线监测会产生大量的误报警或漏报警。

发明内容

为了满足条盒包装机BV在线监测与故障诊断的实际需要，本发明面向生产过程的多工况特性，提出一种包含离线工况识别和划分、稳定工况统计建模、多工况过程在线监测、故障原因诊断等关键技术的超高速条盒包装机在线监测与故障诊断方法。

为了实现上述的目的，本发明采用了以下的技术方案：

卷烟工厂超高速条盒包装机多工况过程在线监测与故障诊断方法，包括：

1)根据超高速条盒包装机BV的运行机理确定统计模型的监测变量，采集某一机台底层PLC中的BV运行数据，获得二维原始数据D(K×J)，K为采样点个数，J为监测变量个数；

2)基于BV车速检测值，采用滑动时间窗口H方法计算稳定度因子SF，当其值较小时(小于等于阈值ζ)表明过程处于稳定工况，反之表明过程处于过渡工况。采用平滑滤波方法对处于稳定工况的滑动时间窗口H内原始数据进行算术平均，获得一组移动平均数据，将移动平均数据中“BV车速”变量值大于阈值ψ的数据作为稳定工况的有效数据D(N×J)，N为有效数据的个数；

3)稳定工况的有效数据被识别出以后，采用自适应k-means聚类方法对有效数据的相似度进行聚类分析，获得每一种稳定工况的建模数据C表示将稳定工况划分成C个聚类，其中第r个聚类为K_rm表示第r个聚类中稳定工况的个数，每个聚类中具有相似过程特性的建模数据可以用一个统计模型来描述；

4)对每一种稳定工况的建模数据(r＝1,2,…,C)进行均值中心化和方差归一化的数据预处理，获得处理后的建模数据X_r(K_rm×J)。对建模数据X_r进行PCA分解，得到对应的主元个数A_r、负载矩阵P_r、控制限和建立对应第r种稳定工况的PCA监测模型；

5)在线监测时，计算当前滑动时间窗口内数据的稳定度因子SF_new，通过与给定阈值ζ的比较确定当前所处的工况类型。若为过渡工况，则将当前时刻和SPE_new两个监测统计量赋值为0；否则采用平滑滤波方法计算当前滑动时间窗口内的移动平均数据，并计算其与各个聚类中心的欧式距离，调用距离最小的聚类所对应的PCA监测模型在线计算和SPE_new两个监测统计量，并与该监测模型的控制限进行比较。如果两个统计量都位于控制限以内，表明设备状态正常，如果至少其中一个统计量超出控制限，表明设备状态异常；

6)当检测到BV有异常时，根据超限统计量处于T²的主元子空间或是SPE的残差子空间，计算每个过程变量对超限统计量的贡献率，贡献率较大的变量被初步确定为造成设备异常的原因变量。

所述步骤2)中，采用滑动时间窗口H方法根据BV车速检测值计算稳定度因子SF，计算公式如下：

其中V_{H max}表示滑动时间窗口H内BV车速的最大值，V_{H min}表示滑动窗口H内BV车速的最小值，V_max表示超高速包装机组允许的车速最大值，H表示滑动时间窗口的长度。

在步骤3)中，采用自适应k-means聚类方法对有效数据的相似度进行聚类分析，通过定义两个聚类中心的最小距离的阈值θ，同时最小化同一聚类数据之间的局部距离平方和以及所有聚类之间的全局距离平方和，获得最优的聚类数目和聚类中心，实现对稳定工况的自适应划分。

在步骤4)中，对建模数据预处理包括依次进行的减均值、除标准差处理获得X_r(K_rm×J)，计算公式如下：

其中下标r代表工况、k代表样本、j代表变量，表示第r种稳定工况第j个变量的均值，表示第r种稳定工况第j个变量的标准差。

将数据X_r进行PCA分解，PCA分解的计算公式为：

其中T_r为主元子空间中(K_rm×A_r)维的得分矩阵；P_r为主元子空间中(J×A_r)维的负载矩阵；E_r为残差子空间中(K_rm×J)维的残差矩阵。A_r表示主元个数，由累计贡献率方法确定。离线计算主元子空间的Hotelling-T²和残差子空间的SPE(Squared PredictionError)监测统计

量。监测统计量的计算公式为：

其中x_rk为建模数据X_r的第k行，P_r为主元子空间的负载矩阵，对角矩阵由建模数据X_r的协方差矩阵∑X_r ^TX_r的前A_r个特征值所构成。

SPE_r监测统计量的计算公式为：

其中表示重构得到的(1×J)维估计向量。

根据F分布计算监测统计量的控制限为：

其中α为置信度，F_α(A_r,K_rm-A_r)是对应于置信度为α，自由度为A_r、K_rm-A_r条件下的F分布临界值。

根据χ²分布计算SPE_r监测统计量的控制限为：

其中g＝v/2n，h＝2n²/v，n、v分别为SPE_r监测统计量的均值和方差。

在步骤5)中，计算处于稳定工况滑动时间窗口的移动平均数据与各个聚类中心的欧式距离，根据距离最小匹配该时刻对应的稳定工况监测模型。欧式距离的计算公式为：

其中x_new,j表示移动平均数据的第j个变量的检测值，Cen_i,j表示第i个聚类中心第j个变量的值，C表示聚类中心的个数。

设当前移动平均数据对应第L(1≤L≤C)个稳定工况监测模型，利用该监测模型建模数据的均值和标准差对移动平均数据进行标准化处理，获得x_Lnew(1×J)。将其投影到第L个监测模型上，计算和SPE_Lnew监测统计量指标。

其中P_L表示第L个监测模型主元子空间的负载矩阵，对角矩阵S_L是由建模数据X_L的协方差矩阵的前A_L个特征值所构成，表示重构得到的(1×J)维估计向量。

将上述监测统计量计算值与第L个监测模型的控制限进行比较，计算公式如下：

如果两个监测统计量都位于控制限内，表明条盒包装机BV处于正常状态，如果其中至少一个监测统计量超出控制限，表明处于异常状态。

在步骤6)中，当任一监测统计量超出控制限时，计算各个过程变量对超限统计量的贡献值，确定引起异常的原因变量。主元子空间监测统计量的超出正常控制限时，主成分得分变量t_La对的贡献率可以计算如下：

其中t_La表示对应第L个监测模型的第a个主成分得分变量，λ_La表示对应第L个监测模型的第a个特征值。

过程变量x_Lnew,j对t_La的贡献率可以计算如下：

其中p_La,j表示对应第L个监测模型的第a个负载向量的第j个变量。

当残差子空间SPE_Lnew监测统计量的超出正常控制限时，过程变量x_Lnew,j对SPE_Lnew的贡献率可以计算如下：

其中表示重构估计向量的第j个变量，表示残差的正负信息。

本发明的有益效果是：

本发明提出一种卷烟工厂超高速条盒包装机多工况过程在线监测与故障诊断方法，主要涉及离线工况识别和划分、稳定工况统计建模、多工况过程在线监测、故障原因辨识等关键技术。基于海量历史数据，通过滑动时间窗口离线计算稳定度因子实现了条盒包装机BV多工况过程稳定工况和过渡工况的有效识别，通过自适应k-means聚类方法对稳定工况自适应划分形成若干个稳定工况数据簇，采用多模型主元分析方法实现每种稳定工况的统计监测建模以及控制限的求取。在线监测时，根据当前滑动时间窗口内数据的稳定度因子判断工况类型，若为过渡工况，则将当前时刻T²、SPE两个监测统计量赋值为0；若为稳定工况，根据稳定工况移动平均数据与各个聚类中心的相似度匹配该时刻对应的监测模型，通过两个监测统计量的在线计算以及与相应控制限的比较，实现超高速条盒包装机BV多工况过程的有效监测。检测到异常发生时，采用变量贡献图方法在统计量超限的子空间中追溯和确定引起故障的原因变量。最终形成一套集工况识别、统计建模、在线监测和故障诊断于一体的切实可行的实施方案，能有效提高超高速条盒包装机BV监测和诊断结果的可靠性和准确性，进而为设备预防维修提供了客观、量化的参考信息。

附图说明

图1为本发明的卷烟工厂超高速条盒包装机多工况过程在线监测与故障诊断流程图。

图2为工况识别前1～1000采样点的BV车速曲线图。

图3为1～1000采样点BV车速的稳定度因子SF曲线图。

图4为经过平滑滤波后稳定工况的BV车速曲线图。

图5为稳定工况有效数据的BV车速曲线图。

图6为BV多工况过程正常状态测试数据的在线监测图。

图7为BV多工况过程故障状态数据1的在线监测图。

图8为BV多工况过程故障状态数据2的在线监测图。

图9为BV多工况过程故障状态数据3的在线监测图。

图10为故障1下主元子空间首次故障报警时刻的T²变量贡献图。

图11为故障1下残差子空间首次故障报警时刻的SPE变量贡献图。

图12为故障2下残差子空间首次故障报警时刻的SPE变量贡献图。

图13为故障3下主元子空间首次故障报警时刻的T²变量贡献图。

具体实施方式

为了更好的理解本发明的技术方案，以下结合说明书附图对本发明的实施方式作进一步描述。卷烟包装过程中G.D公司超高速包装设备由小盒包装机H1000、小盒薄膜封装机W1000和条盒包装机BV三部分组成，本实施是用于超高速条盒包装机BV的在线监测与故障诊断方法。条盒包装机BV是一台高速自动外条包装机，可以包装各种不同特性和硬度的条包。本发明的卷烟工厂超高速条盒包装机多工况过程在线监测与故障诊断方法的实现框图如图1所示，方法主要分为以下几步：

(1)历史运行数据的获取

根据超高速条盒包装机BV的运行机理确定统计模型的监测变量，采集某一机台底层PLC中的BV运行数据，获得二维原始数据D(K×J)，K为采样点个数，J为监测变量个数。

本实例中，选取A4机台的BV条盒包装机，采集底层PLC中的BV运行数据，包含以下变量：伺服电机电流、伺服电机温度、加热器温度、BV车速等43个监测变量，如表1所示。监测变量的采样频率为10秒/次，获得的原始数据为D(56322×43)。

表1超高速条盒包装机BV的建模变量表

(2)离线工况的识别

基于BV车速检测值，采用滑动时间窗口H方法计算稳定度因子SF，当其值较小时(小于等于阈值ζ)表明过程处于稳定工况，反之表明过程处于过渡工况。采用平滑滤波方法对处于稳定工况的滑动时间窗口H内原始数据进行算术平均，获得一组移动平均数据，将移动平均数据中“BV车速”变量值大于阈值ψ的数据作为稳定工况的有效数据D(N×J)，N为有效数据的个数。

采用滑动时间窗口H方法计算稳定度因子SF，计算公式如下：

其中V_Hmax表示滑动时间窗口H内BV车速的最大值，V_Hmin表示滑动窗口H内BV车速的最小值，V_max表示超高速包装机组允许的车速最大值，H表示滑动时间窗口的长度。

本实例中，滑动时间窗口的长度H取为12，允许的车速最大值V_max取为1000，阈值ζ取为0.01，阈值ψ取为80，获得的有效数据为D(22260×43)。以前1000个原始采样数据为例，工况识别前BV车速如图2所示，相应的稳定度因子如图3所示，经过平滑滤波后稳定工况的BV车速如图4所示，稳定工况有效数据的BV车速如图5所示。

(3)离线工况的划分

稳定工况的有效数据被识别出以后，采用自适应k-means聚类方法对有效数据的相似度进行聚类分析，获得每一种稳定工况的建模数据C表示将稳定工况划分成C个聚类，每个聚类中具有相似过程特性的建模数据可以用一个统计模型来描述。自适应k-means聚类方法的求取过程如下：

步骤1：设定初始聚类中心个数C₀，两个聚类中心的最小距离阈值θ，每个聚类俘获的最小数据量阈值minL。初始化C₀个聚类中心W_i,c(i＝1,…,J；c＝1,…,C₀)；

步骤2：计算两个聚类中心的欧式距离dist(W_i,c1，W_i,c2)，若其小于预定的阈值θ，剔除其中一个聚类中心；

步骤3：计算每个稳定工况的有效数据D_k(k＝1,2,…,N)到所有聚类中心的欧式距离dist(D_k,W_i,c)，若D_k与第c^*类的中心的距离最小，则将D_k的隶属关系定义为m(k)＝c^*；

步骤4：I_num次迭代后，若某聚类中心俘获的数据样本个数小于预定的阈值minL，则剔除该类。

步骤5：更新聚类个数为C_i+1，并根据数据的隶属关系重新计算新的聚类中心W_i+1,c(c＝1,2,…,C_i+1)。

步骤6：如果算法满足收敛条件则结束，否则返回步骤2，进行下一次迭代计算。收敛条件有：两次迭代中的聚类中心距离的变化小于一个很小的阈值ε，或者每个聚类中数据D_k到聚类中心的距离平方和以及聚类之间的距离平方和达到最小。

本实例中，初始聚类中心个数C₀取为20，最小距离阈值θ取为150，最小数据量阈值minL取为500，阈值ε取为0.0001。通过自适应聚类，稳定工况有效数据被划分为9个聚类(即C＝9)，获得每一种稳定工况的建模数据

(4)稳定工况的PCA统计建模

对每一种稳定工况的建模数据(r＝1，2，…，C)进行均值中心化和方差归一化的数据预处理，获得处理后的建模数据X_r(K_rm×J)。对建模数据X_r进行PCA分解，得到对应的主元个数A_r、负载矩阵P_r、控制限和建立对应第r个稳定工况的PCA监测模型。

对建模数据预处理包括依次进行的减均值、除标准差处理获得X_r(K_rm×J)，计算公式如下：

其中下标r代表工况数、k代表样本数、j代表变量数，表示第r种稳定工况第j个变量的均值，表示第r种稳定工况第j个变量的标准差。

将数据X_r进行PCA分解，PCA分解的计算公式为：

其中T_r为主元子空间中(K_rm×A_r)维的得分矩阵；P_r为主元子空间中(J×A_r)维的负载矩阵；E_r为残差子空间中(K_rm×J)维的残差矩阵。A_r表示主元个数，由累计贡献率方法确定。

离线计算主元子空间的Hotelling-T²和残差子空间的SPE(Squared PredictionError)监测统计量。监测统计量的计算公式为：

其中x_rk为建模数据X_r的第k行，P_r为主元子空间的负载矩阵，对角矩阵由建模数据X_r的协方差矩阵ΣX_r ^TX_r的前A_r个特征值所构成。

SPE_r监测统计量的计算公式为：

其中表示重构得到的(1×J)维估计向量。

根据F分布计算监测统计量的控制限为：

其中α为置信度，F_α(A_r，K_rm-A_r)是对应于置信度为α，自由度为A_r、K_rm-A_r条件下的F分布临界值。

根据χ²分布计算SPE_r监测统计量的控制限为：

其中g＝v/2n，h＝2n²/v，n、v分别为SPE_r监测统计量的均值和方差。

本实例中，9个PCA监测模型的主元个数依次为A₁＝27、A₂＝28、A₃＝28、A₄＝28、A₅＝27、A₆＝28、A₇＝27、A₈＝28、A₉＝27，则P₁为(43×27)负载矩阵、P₂为(43×28)负载矩阵、P₃为(43×28)负载矩阵、P₄为(43×28)负载矩阵、P₅为(43×27)负载矩阵、P₆为(43×28)负载矩阵、P₇为(43×27)负载矩阵、P₈为(43×28)负载矩阵、P₉为(43×27)负载矩阵。统计控制限的置信度α取为0.99，监测模型1的控制限和监测模型2的控制限和监测模型3的控制限和监测模型4的控制限和监测模型5的控制限和ctr_SPE5＝13.75、监测模型6的控制限和监测模型7的控制限和监测模型8的控制限和监测模型9的控制限和

(5)多工况过程的在线监测

在线监测时，计算当前滑动时间窗口内数据的稳定度因子SF_new，通过与给定阈值ζ的比较确定当前所处的工况类型。若为过渡工况，则将当前时刻和SPE_new两个监测统计量赋值为0；否则采用平滑滤波方法计算当前滑动时间窗口内的移动平均数据，并计算其与各聚类中心的欧式距离，调用距离最小的聚类所对应的PCA监测模型在线计算和SPE_new两个监测统计量。如果两个统计量都位于该监测模型对应的控制限以内，表明设备状态正常，如果至少其中一个统计量超出控制限，表明设备状态异常。

计算处于稳定工况滑动时间窗口的移动平均数据与各个聚类中心的欧式距离，根据距离最小匹配该时刻对应的稳定工况监测模型。欧式距离的计算公式为：

其中x_new,j表示移动平均数据的第j个变量的检测值，Cen_i,j表示第i个聚类中心第j个变量的值，C表示聚类中心的个数。

设当前移动平均数据对应第L(1≤L≤C)个稳定工况监测模型，利用该监测模型建模数据的均值和标准差对移动平均数据进行标准化处理，获得x_Lnew(1×J)。将其投影到第L个监测模型上，计算和SPE_Lnew监测统计量指标。

其中P_L表示第L个监测模型主元子空间的负载矩阵，对角矩阵S_L是由建模数据X_L的协方差矩阵的前A_L个特征值所构成，表示重构得到的(1×J)维估计向量。

将上述监测统计量计算值与第L个监测模型的控制限进行比较，计算公式如下：

如果两个监测统计量都位于控制限内，表明条盒包装机BV处于正常状态，如果其中至少一个监测统计量超出控制限，表明处于异常状态。

本实例中，选取条盒包装机BV多工况过程中的8941组采样数据D _new(8941×43)进行分析，该测试数据T²和SPE统计量的监测结果如图6所示。引入首次故障报警时间，定义为稳定工况下，连续9个监测统计量超过对应控制限的采样时间点。在两个监测子空间中有个别样本点的监测统计量超出控制限，根据首次故障报警的定义，主元子空间中有3个故障报警、残差子空间中未发生故障报警，基本能够对多工况过程的正常状态做出正确指示，表明所提在线监测与故障诊断方法具有准确监测多工况过程正常状态的能力。

作为比较，重新选取条盒包装机BV多工况过程中的9687组采样数据D _new(9687×43)进行异常状态在线监测与故障诊断的分析。在第4501个样本点的时候，分别令A914B电机电流阶跃增长25％(突变故障)、下方加热器二温度阶跃降低5％(突变故障)、侧边内侧再次加热器温度以0.002的斜率增加(渐变故障)，一直持续到第6500个样本点的时候恢复正常，形成故障状态数据1、2和3。故障状态数据T²和SPE统计量的监测结果分别如图7至图9所示。

对于故障1，主元子空间的T²统计量和残差子空间的SPE统计量基本同时出现明显的连续超限，主元子空间、残差子空间检测出的首次故障报警时间分别为第4533、4532个样本点时刻，即监测模型延后310秒准确检测出故障的发生；对于故障2，残差子空间的SPE统计量出现明显的连续超限，检测出的首次故障报警时间为第4532个样本点时刻，即监测模型延后310秒准确检测出故障的发生；对于故障3，主元子空间的T²统计量先于残差子空间的SPE统计量出现明显的连续超限，主元子空间检测出的首次故障报警时间为第4898个样本点时刻，即监测模型延后3970秒准确检测出故障的发生。说明所提方法能对多工况过程的突变和渐变故障做出正确指示，表明该方法具有准确监测多工况过程故障状态的能力。

(6)基于贡献图的故障诊断

当任一监测统计量超出控制限时，计算各个过程变量对超限统计量的贡献值，确定引起异常的原因变量。主元子空间监测统计量的超出正常控制限时，主成分得分变量t_La对的贡献率可以计算如下：

其中t_La表示对应第L个监测模型的第a个主成分得分变量，λ_La表示对应第L个监测模型的第a个特征值。

过程变量x_Lnew，j对t_La的贡献率可以计算如下：

其中p_La,j表示对应第L个监测模型的第a个负载向量的第j个变量。

当残差子空间SPE_Lnew监测统计量的超出正常控制限时，过程变量x_Lnew,j对SPE_Lnew的贡献率可以计算如下：

其中表示重构估计向量的第j个变量，表示残差的正负信息。

本实例中，根据图7的监测结果，故障发生在主元子空间和残差子空间。在主元子空间首次故障报警的第4533样本点时刻，过程变量对监测统计量超限的贡献率如图10所示，在残差子空间首次故障报警的第4532个样本点时刻，过程变量对SPE_Lnew监测统计量超限的贡献率如图11所示。确定引起故障的原因变量均为第7个变量，即A914B电机电流。

根据图8的监测结果，故障发生在残差子空间。在残差子空间首次故障报警的第4532个样本点时刻，过程变量对SPE_Lnew监测统计量超限的贡献率如图12所示。确定引起故障的原因变量为第34个变量，即下方加热器二温度。

根据图9的监测结果，故障发生在主元子空间和残差子空间。在较早检测出故障的主元子空间中，首次故障报警的第4898样本点时刻，过程变量对监测统计量超限的贡献率如图13所示。确定引起故障的原因变量为第40个变量，即侧边内侧再次加热器温度。

Claims (6)

1.卷烟工厂超高速条盒包装机多工况过程在线监测与故障诊断方法，其特征在于该方法包括以下的步骤：

1)根据超高速条盒包装机BV的运行机理确定统计模型的监测变量，采集某一机台底层PLC中的BV运行数据，获得二维原始数据D(K×J)，K为采样点个数，J为监测变量个数；

2)基于BV车速检测值，采用滑动时间窗口H方法计算稳定度因子SF，当其值小于等于阈值ζ表明过程处于稳定工况，反之表明过程处于过渡工况；采用平滑滤波方法对处于稳定工况的滑动时间窗口H内原始数据进行算术平均，获得一组移动平均数据，将移动平均数据中“BV车速”变量值大于阈值ψ的数据作为稳定工况的有效数据D(N×J)，N为有效数据的个数；

3)稳定工况的有效数据被识别出以后，采用自适应k-means聚类方法对有效数据的相似度进行聚类分析，获得每一种稳定工况的建模数据C表示将稳定工况划分成C个聚类，其中第r个聚类为K_rm表示第r个聚类中稳定工况的个数，每个聚类中具有相似过程特性的建模数据可以用一个统计模型来描述；

4)对每一种稳定工况的建模数据进行均值中心化和方差归一化的数据预处理，获得处理后的建模数据X_r(K_rm×J)；对建模数据X_r进行PCA分解，得到对应的主元个数A_r、负载矩阵P_r、控制限和建立对应第r种稳定工况的PCA监测模型；

5)在线监测时，计算当前滑动时间窗口内数据的稳定度因子SF_new，通过与给定阈值ζ的比较确定当前所处的工况类型；若为过渡工况，则将当前时刻和SPE_new两个监测统计量赋值为0；否则采用平滑滤波方法计算当前滑动时间窗口内的移动平均数据，并计算其与各个聚类中心的欧式距离，调用距离最小的聚类所对应的PCA监测模型在线计算和SPE_new两个监测统计量，并与该监测模型的控制限进行比较；如果两个统计量都位于控制限以内，表明设备状态正常，如果至少其中一个统计量超出控制限，表明设备状态异常；

6)当检测到BV有异常时，根据超限统计量处于T²的主元子空间或是SPE的残差子空间，计算每个过程变量对超限统计量的贡献率，贡献率较大的变量被初步确定为造成设备异常的原因变量。

2.根据权利要求1所述的卷烟工厂超高速条盒包装机多工况过程在线监测与故障诊断方法，其特征在于所述步骤2)中，采用滑动时间窗口H方法根据BV车速检测值计算稳定度因子SF，计算公式如下：

其中V_{H max}表示滑动时间窗口H内BV车速的最大值，V_{H min}表示滑动窗口H内BV车速的最小值，V_max表示超高速包装机组允许的车速最大值，H表示滑动时间窗口的长度。

3.根据权利要求1所述的卷烟工厂超高速条盒包装机多工况过程在线监测与故障诊断方法，其特征在于在步骤3)中，采用自适应k-means聚类方法对有效数据的相似度进行聚类分析，通过定义两个聚类中心的最小距离的阈值θ，同时最小化同一聚类数据之间的局部距离平方和以及所有聚类之间的全局距离平方和，获得最优的聚类数目和聚类中心，实现对稳定工况的自适应划分。

4.根据权利要求1所述的卷烟工厂超高速条盒包装机多工况过程在线监测与故障诊断方法，其特征在于在步骤4)中，对建模数据预处理包括依次进行的减均值、除标准差处理获得X_r(K_rm×J)，计算公式如下：

其中下标r代表工况、k代表样本、j代表变量，表示第r种稳定工况第j个变量的均值，表示第r种稳定工况第j个变量的标准差；

将数据X_r进行PCA分解，PCA分解的计算公式为：

其中T_r为主元子空间中(K_rm×A_r)维的得分矩阵；P_r为主元子空间中(J×A_r)维的负载矩阵；E_r为残差子空间中(K_rm×J)维的残差矩阵；A_r表示主元个数，由累计贡献率方法确定；

离线计算主元子空间的Hotelling-T2和残差子空间的SPE(Squared PredictionError)监测统计量；监测统计量的计算公式为：

其中x_rk为建模数据X_r的第k行，P_r为主元子空间的负载矩阵，对角矩阵由建模数据X_r的协方差矩阵∑X_r ^TX_r的前A_r个特征值所构成；

SPE_r监测统计量的计算公式为：

其中表示重构得到的(1×J)维估计向量；

根据F分布计算监测统计量的控制限为：

其中α为置信度，F_α(A_r,K_rm-A_r)是对应于置信度为α，自由度为A_r、K_rm-A_r条件下的F分布临界值；

根据χ²分布计算SPE_r监测统计量的控制限为：

其中g＝v/2n，h＝2n²/v，n、v分别为SPE_r监测统计量的均值和方差。

5.根据权利要求1所述的卷烟工厂超高速条盒包装机多工况过程在线监测与故障诊断方法，其特征在于在步骤5)中，计算处于稳定工况滑动时间窗口的移动平均数据与各个聚类中心的欧式距离，根据距离最小匹配该时刻对应的稳定工况监测模型；欧式距离的计算公式为：

其中x_new,j表示移动平均数据的第j个变量的检测值，Cen_i,j表示第i个聚类中心第j个变量的值，C表示聚类中心的个数；

设当前移动平均数据对应第L个稳定工况监测模型，1≤L≤C，利用该监测模型建模数据的均值和标准差对移动平均数据进行标准化处理，获得x_Lnew(1×J)；将其投影到第L个监测模型上，计算和SPE_Lnew监测统计量指标；

其中P_L表示第L个监测模型主元子空间的负载矩阵，对角矩阵S_L是由建模数据X_L的协方差矩阵的前A_L个特征值所构成，表示重构得到的(1×J)维估计向量；

将上述监测统计量计算值与第L个监测模型的控制限进行比较，计算公式如下：

如果两个监测统计量都位于控制限内，表明条盒包装机BV处于正常状态，如果其中至少一个监测统计量超出控制限，表明处于异常状态。

6.根据权利要求1所述的卷烟工厂超高速条盒包装机多工况过程在线监测与故障诊断方法，其特征在于在步骤6)中，当任一监测统计量超出控制限时，计算各个过程变量对超限统计量的贡献值，确定引起异常的原因变量；主元子空间监测统计量的超出正常控制限时，主成分得分变量t_La对的贡献率可以计算如下：

其中t_La表示对应第L个监测模型的第a个主成分得分变量，λ_La表示对应第L个监测模型的第a个特征值；

过程变量x_Lnew,j对t_La的贡献率可以计算如下：

其中p_La,j表示对应第L个监测模型的第a个负载向量的第j个变量；

当残差子空间SPE_Lnew监测统计量的超出正常控制限时，过程变量x_Lnew,j对SPE_Lnew的贡献率可以计算如下：

其中表示重构估计向量的第j个变量，表示残差的正负信息。