""" @Time : 2024-03-04 @Author : Wei @FileName: DiagnosisLib.py @Software: 故障诊断函数库 """ import numpy as np import matplotlib.pyplot as plt from math import pi from scipy import signal from scipy import interpolate import math # git外网使用 class DiagnosisLib: def __init__(self): pass ############################################################################### ############################################################################### # 半带滤波器系数 def halfbandcoefficient(self): halfbandcoefficient = np.array([0.00000063557331607068,0,-0.00000627255234261115,0,0.00002474955105003003,0,\ -0.00007140981582044673,0,0.00017237915446870680,0,-0.00036863887543283979,0,\ 0.00072039609678193576,0,-0.00131147951710543810,0,0.00225360507770744260,0,\ -0.00369069002024503300,0,0.00580407535380912750,0,-0.00882079299392407940,0,\ 0.01302944110367636400,0,-0.01881320502983020000,0,0.02672090662754737100,0,\ -0.03762627155576255700,0,0.05311294860790129200,0,-0.07653440360979835200,0,\ 0.11663131038150232000,0,-0.20562506738991029000,0,0.63439761518650584000,1,\ 0.63439761518650584000,0,-0.20562506738991029000,0,0.11663131038150232000,0,\ -0.07653440360979835200,0,0.05311294860790129200,0,-0.03762627155576255700,0,\ 0.02672090662754737100,0,-0.01881320502983020000,0,0.01302944110367636400,0,\ -0.00882079299392407940,0,0.00580407535380912750,0,-0.00369069002024503300,0,\ 0.00225360507770744260,0,-0.00131147951710543810,0,0.00072039609678193576,0,\ -0.00036863887543283979,0,0.00017237915446870680,0,-0.00007140981582044673,0,\ 0.00002474955105003003,0,-0.00000627255234261115,0,0.00000063557331607068]) return halfbandcoefficient ############################################################################### ############################################################################### # 计算时域波形均值 def Average(self,xTimeWave): return sum(xTimeWave) / len(xTimeWave) ############################################################################### ############################################################################### # 计算时域波形xTimeWave的频谱 # fs是时域波形的采样频率，单位是Hz # 调用代码须保证输入时域波形xTimeWave的数据点数为2的N幂 # 频谱计算默认采用汉宁窗 # 输出为频谱的x轴freq即频率，单位是Hz， 和频谱的y轴，单位是时域波形的幅值单位 # x轴和y轴的数据点数为输入时域波形xTimeWave数据点数的1/2.56倍+1，即1024点时域波形 # 输出频谱的x轴和y轴分别有401点数据，2048点时域波形时，输出频谱的x轴和y轴分别有801点数据。 def Spectrum(self,xTimeWave,fs): DataLen = len(xTimeWave) xSpec = abs(np.fft.fft((xTimeWave*np.hanning(DataLen))))/DataLen*4 freq = np.arange(0, fs / 2.56 + fs / DataLen, fs / DataLen) return [freq, xSpec[0:int(DataLen/2.56)+1]] ############################################################################### ############################################################################### # 计算时域波形xTimeWave的幅值谱倒频谱 # fs是时域波形的采样频率，单位是Hz # 调用代码须保证输入时域波形xTimeWave的数据点数为2的N幂 # 输出为倒谱的x轴quefrency,单位是秒，倒谱的幅值单位是dB def Cepstrum(self,xTimeWave,fs): # cepstrum = np.abs(np.fft.ifft(np.log(np.absolute(np.fft.fft(xTimeWave))))) cepstrum = 20*np.log((np.abs(np.fft.ifft(np.log(np.abs(np.fft.fft(xTimeWave))))))/0.1) quefrency = np.array(range(len(xTimeWave)))/fs return [quefrency,cepstrum] ############################################################################### ############################################################################### # 计算频谱xSpec的固定频率特征值,xSpec是频谱的y轴数组，其数据点数为101，201，401，801，1601，3201，6401，12801 # fMax是频谱的分析频宽，单位为Hz # fTarget是特征频率的目标值，单位为Hz。比如轴的转频，齿轮啮合频率，水泵叶片通过频率等等 # fRange是特征频率的搜索范围，其输入值根据fRangeMode的不同值而不一样， # fRangeMode可以赋值0或1，默认值是0。当fRangeMode = 0时，fRange的值是xx%， 比如 2%(实际输入为0.02) ， # 此时搜索范围是 fTarge +/- fTarge*fRange/2。 # 当fRangeMode = 1时，fRange的值是频率值，比如2Hz， # 此时的搜索范围就是 fTarge +/- fRange/2 # Detection是特征值的峰值（此时输入0），有效值（输入1），峰峰值（输入2） # numHarmonics为谐波数，0代表不需要考虑谐波即只有特征频率的基频，1代表包括基频和2倍频，2代表包括基频和2、3倍频 def FixedFreqFeature(self,xSpec, fMax, fTarget, fRange, fRangeMode, Detection, numHarmonics): if fRangeMode == '百分比': fRange = fRange*fTarget if fRangeMode == '数值': fRange = fRange # print(fRange) SpecResolution = fMax/(len(xSpec)-1) index_SearchRangeMax = int(math.ceil((fTarget + fRange/2)/SpecResolution)) index_SearchRangeMin = int(math.floor((fTarget - fRange/2)/SpecResolution)) if index_SearchRangeMin < 3: index_SearchRangeMin = 3 index_SearchRange = index_SearchRangeMax - index_SearchRangeMin xArray = xSpec[index_SearchRangeMin:(index_SearchRangeMax+1)] while((int(np.argmax(xArray))==0) or (int(np.argmax(xArray))==len(xArray)-1)): index_SearchRangeMin = index_SearchRangeMin -1 index_SearchRange = index_SearchRange + 1 xArray = xSpec[int(index_SearchRangeMin-math.floor(index_SearchRange/2)):int(index_SearchRangeMin+math.ceil(index_SearchRange/2+1))] ActualTargetFreq = self.FundamentalFreqActualValue(xArray, SpecResolution,index_SearchRangeMin) # print(ActualTargetFreq) featureValue = self.FreqFeatureValue(xSpec, SpecResolution, ActualTargetFreq, index_SearchRange, Detection, numHarmonics) return featureValue ############################################################################### ############################################################################### # 该函数只用于内部调用，其调用函数为 FixedFreqFeature def FundamentalFreqActualValue(self,xArray, SpecResolution, index_SearchRangeMin): indexPeak = int(self.PeakFinderInArray(xArray)) indexCenter = int(indexPeak + index_SearchRangeMin) ActualFundamentalFreq = self.FrequencyCompensation(indexCenter,xArray[indexPeak],xArray[indexPeak-1],xArray[indexPeak+1],SpecResolution) return ActualFundamentalFreq ############################################################################### ############################################################################### # 该函数只用于内部调用，其调用函数为 FundamentalFreqActualValue def PeakFinderInArray(self,xArray): NoElements = len(xArray) index_peak = int(np.argmin(xArray)) # initialize the peak index with the minimal value's index peak = xArray[index_peak] # initialize the peak value with the minimal value for i in range(0, NoElements-2): if (xArray[i+1] > xArray[i]) and (xArray[i+1] > xArray[i+2]) and (xArray[i+1] > peak): index_peak = i+1 peak = xArray[i+1] return index_peak ############################################################################### ############################################################################### # 该函数只用于内部调用，其调用函数为 FundamentalFreqActualValue def FrequencyCompensation(self,indexCenter, AmpCenter, AmpLeft, AmpRight, xSpecResolution): if AmpRight >= AmpLeft: delta_k = (2*AmpRight-AmpCenter)/(AmpRight+AmpCenter) else: delta_k = (AmpCenter-2*AmpLeft)/(AmpLeft+AmpCenter) ActualFrequency = (indexCenter + delta_k)*xSpecResolution return ActualFrequency ############################################################################### ############################################################################### # 该函数只用于内部调用，其调用函数为 FixedFreqFeature def FreqFeatureValue(self,xSpec, SpecResolution, ActualTargetFreq, index_SearchRange, Detection, numHarmonics = 0): sumation = 0 # print(ActualTargetFreq) featureValue = [] featureValue = [0 for i in range(numHarmonics+1)] for i in range(0,numHarmonics+1): index_i = int(round((i+1)*ActualTargetFreq/SpecResolution)) a = xSpec[int(index_i-math.floor(index_SearchRange/2)):int(index_i+math.ceil(index_SearchRange/2+1))] while((int(np.argmax(a))==0) or (int(np.argmax(a))==len(a)-1)): index_SearchRange = index_SearchRange + 2 a = xSpec[int(index_i-math.floor(index_SearchRange/2)):int(index_i+math.ceil(index_SearchRange/2+1))] indexPeak = int(self.PeakFinderInArray(a) + index_i-math.floor(index_SearchRange/2)) # print(int(index_i-math.floor(index_SearchRange/2)),indexPeak,int(index_i+math.ceil(index_SearchRange/2))) for j in range(indexPeak-2,indexPeak+2+1): sumation = sumation + xSpec[j]*xSpec[j] featureValue[i] = np.sqrt(sumation/1.5) sumation = 0 # print(featureValue[i]) if Detection == '峰值': pass if Detection == '有效值': featureValue = [x/math.sqrt(2) for x in featureValue] if Detection == '峰峰值': featureValue = [x*2 for x in featureValue] return featureValue ############################################################################### ############################################################################### # 计算频域积分 # xSpecOrig需要对其计算积分的频谱数组，其长度为：101，201，401，801，1601，3201，6401，12801 # fMax为频谱的分析频宽，单位为Hz # fCut为计算积分的低频截止频率，单位为Hz def FreqIntegration(self,xSpecOrig, fMax, fCut): nLine = len(xSpecOrig) - 1 lineResolution = fMax / nLine iStarting = int(math.floor(fCut / lineResolution)) if iStarting == 0: iStarting = 1 xSpecOut = np.zeros(len(xSpecOrig)) for i in range(iStarting, len(xSpecOrig)): xSpecOut[i] = xSpecOrig[i]/(2*pi*lineResolution*i) return xSpecOut ############################################################################### ############################################################################### # 巴特沃斯高通滤波 # xTimeWave是时域波形数组， fs是采样频率，单位为Hz # fco 为高通截止频率，单位为Hz # Order为滤波器阶数，通常为偶数，2，4，6，8一般不超过8 def HighPassButter(self,xTimeWave, fs, fco, Order): sos = signal.butter(Order, fco/(fs/2), 'hp', output='sos') filtered = signal.sosfilt(sos, xTimeWave) return filtered ############################################################################### ############################################################################### # 巴特沃斯低通滤波 # xTimeWave是时域波形数组， fs是采样频率，单位为Hz # fco 为低通截止频率，单位为Hz # Order为滤波器阶数，通常为偶数，2，4，6，8一般不超过8 def LowPassButter(self,xTimeWave, fs, fco, Order): sos = signal.butter(Order, fco/(fs/2), 'lp', output='sos') filtered = signal.sosfilt(sos, xTimeWave) return filtered ############################################################################### ############################################################################### # 计算包络 # xTimeWave是时域波形数组 # DownSampleRate是降采样率，调用函数需提供此数值，其计算公式为 # 原始波形xTimeWave的分析频宽/包络后时域波形的分析频宽，比如原始波形分析频宽是10000Hz(采样频率是25600Hz),包络后时域 # 波形的分析频宽是1000Hz，那 DownSampleRate = 10000/1000 = 10 # Method是SKF的包络方法，还是PeakHold包络方法 def EnvLopInTime(self,xTimeWave, DownSampleRate, Method = "SKF"): xOut = np.zeros(math.floor(len(xTimeWave)/DownSampleRate)) if Method == "SKF": scale = 2.2 xTimeWave = np.absolute(scale*xTimeWave) w = 1/(1.28*DownSampleRate) Order = 4 sos = signal.butter(Order, w, 'lp', output='sos') xTimeWave = signal.sosfilt(sos, xTimeWave) xOut =xTimeWave[::DownSampleRate] if Method == "PeakHold": for i in range(0, len(xOut)): xOut[i] = np.amax(np.absolute(xTimeWave[i*DownSampleRate:((i+1)*DownSampleRate-1)])) return xOut - np.mean(xOut) ############################################################################### ############################################################################### # 内部调用函数 def DownSampleBy2(self,xTimeWave): xTimeWave = signal.resample_poly(xTimeWave, 1, 2) return xTimeWave def DownSampleBy4(self,xTimeWave): xTimeWave = signal.resample_poly(xTimeWave,1,4) return xTimeWave def DownSampleBy5(self,xTimeWave): xTimeWave = signal.resample_poly(xTimeWave,1,5) return xTimeWave def DownSampleBy8(self,xTimeWave): xTimeWave = signal.resample_poly(xTimeWave,1,8) return xTimeWave def DownSampleBy10(self,xTimeWave): xTimeWave = signal.resample_poly(xTimeWave,1,10) return xTimeWave def DownSampleBy20(self,xTimeWave): xTimeWave = signal.resample_poly(xTimeWave,1,10) xTimeWave = signal.resample_poly(xTimeWave,1,2) return xTimeWave def DownSampleBy40(self,xTimeWave): xTimeWave = signal.resample_poly(xTimeWave,1,10) xTimeWave = signal.resample_poly(xTimeWave,1,4) return xTimeWave def DownSampleBy50(self,xTimeWave): xTimeWave = signal.resample_poly(xTimeWave,1,10) xTimeWave = signal.resample_poly(xTimeWave,1,5) return xTimeWave def DownSampleBy80(self,xTimeWave): xTimeWave = signal.resample_poly(xTimeWave,1,10) xTimeWave = signal.resample_poly(xTimeWave,1,8) return xTimeWave def DownSampleBy100(self,xTimeWave): xTimeWave = signal.resample_poly(xTimeWave,1,10) xTimeWave = signal.resample_poly(xTimeWave,1,10) return xTimeWave def DownSampleBy200(self,xTimeWave): xTimeWave = signal.resample_poly(xTimeWave,1,10) xTimeWave = signal.resample_poly(xTimeWave,1,10) xTimeWave = signal.resample_poly(xTimeWave,1,2) return xTimeWave def DownSampleBy400(self,xTimeWave): xTimeWave = signal.resample_poly(xTimeWave,1,10) xTimeWave = signal.resample_poly(xTimeWave,1,10) xTimeWave = signal.resample_poly(xTimeWave,1,4) return xTimeWave def DownSampleBy500(self,xTimeWave): xTimeWave = signal.resample_poly(xTimeWave,1,10) xTimeWave = signal.resample_poly(xTimeWave,1,10) xTimeWave = signal.resample_poly(xTimeWave,1,5) return xTimeWave def DownSampleBy800(self,xTimeWave): xTimeWave = signal.resample_poly(xTimeWave,1,10) xTimeWave = signal.resample_poly(xTimeWave,1,10) xTimeWave = signal.resample_poly(xTimeWave,1,8) return xTimeWave def DownSampleBy1000(self,xTimeWave): xTimeWave = signal.resample_poly(xTimeWave,1,10) xTimeWave = signal.resample_poly(xTimeWave,1,10) xTimeWave = signal.resample_poly(xTimeWave,1,10) return xTimeWave def DownSampleBy2000(self,xTimeWave): xTimeWave = signal.resample_poly(xTimeWave,1,10) xTimeWave = signal.resample_poly(xTimeWave,1,10) xTimeWave = signal.resample_poly(xTimeWave,1,10) xTimeWave = signal.resample_poly(xTimeWave,1,2) return xTimeWave def DownSampleBy4000(self,xTimeWave): xTimeWave = signal.resample_poly(xTimeWave,1,10) xTimeWave = signal.resample_poly(xTimeWave,1,10) xTimeWave = signal.resample_poly(xTimeWave,1,10) xTimeWave = signal.resample_poly(xTimeWave,1,4) return xTimeWave def DownSample(self, xTimeWave, fsIn, fsOut, nPointIn, nPointOut): #输入: xTimeWave: 时域波形数据 # Fs: xTimeWave波形数据的采样频率(Hz) # fs: 降采样后时域波形数据的采样频率(Hz) # nPointIn: 输入波形的数据点数 # nPointOut: 返回的波形的数据点数 # #返回： 返回一个长度为nPointOut的时域波形数组 DownSampleRate = int(fsIn/fsOut) if int(nPointIn/nPointOut) < DownSampleRate: print('Input data length is not enough') return -1 if DownSampleRate == 2: return self.DownSampleBy2(xTimeWave[0:int(nPointOut*2)]) if DownSampleRate == 4: return self.DownSampleBy4(xTimeWave[0:int(nPointOut*4)]) if DownSampleRate == 5: return self.DownSampleBy5(xTimeWave[0:int(nPointOut*5)]) if DownSampleRate == 8: return self.DownSampleBy8(xTimeWave[0:int(nPointOut*8)]) if DownSampleRate == 10: return self.DownSampleBy10(xTimeWave[0:int(nPointOut*10)]) if DownSampleRate == 20: return self.DownSampleBy20(xTimeWave[0:int(nPointOut*20)]) if DownSampleRate == 40: return self.DownSampleBy40(xTimeWave[0:int(nPointOut*40)]) if DownSampleRate == 50: return self.DownSampleBy50(xTimeWave[0:int(nPointOut*50)]) if DownSampleRate == 100: return self.DownSampleBy100(xTimeWave[0:int(nPointOut*100)]) if DownSampleRate == 200: return self.DownSampleBy200(xTimeWave[0:int(nPointOut*200)]) if DownSampleRate == 400: return self.DownSampleBy400(xTimeWave[0:int(nPointOut*400)]) if DownSampleRate == 500: return self.DownSampleBy500(xTimeWave[0:int(nPointOut*500)]) if DownSampleRate == 1000: return self.DownSampleBy1000(xTimeWave[0:int(nPointOut*1000)]) if DownSampleRate == 2000: return self.DownSampleBy2000(xTimeWave[0:int(nPointOut*2000)]) if DownSampleRate == 4000: return self.DownSampleBy4000(xTimeWave[0:int(nPointOut*4000)]) ############################################################################### ############################################################################### def OrderAnalysis(self, xTimeWave, Fs, xTimeSpeedPulse, nPulsePerRev, nRev, nLine, TransLevel, TransLevelPara): #输入: xTimeWave: 时域波形数据 # Fs: xTimeWave波形数据的采样频率(Hz) # xTimeSpeedPulse: 相邻转速脉冲信号的时间间隔（s） # nPulsePerRev: 每周的脉冲个数 # nRev: 阶次分析整周数（整周数必须是2的N次幂，比如：4、8、16、32、64、128、256、512等等） # nLine： 阶次频谱的谱线数（谱线数可选项为：100、200、400、800、1600、3200、6400） # TransLevel： 阶次分析轴与转速测量轴的转速等级关系， # 0：表示阶次分析轴与转速分析轴是同一根轴 # 1：表示阶次分析轴与转速分析轴不在同一根轴上，有1次转速比变化 # 2：表示阶次分析轴与转速分析轴不在同一根轴上，有2次转速比变化 # 3：表示阶次分析轴与转速分析轴不在同一根轴上，有3次转速比变化， 以此类推 # TransLevelPara： 齿轮齿数数组， # 当 TransLevel == 0 时，这个数组只有1个元素，可以是[0] # 当 TransLevel == 1 时，有一对齿轮啮合，2根轴。这个数组有2个元素，第一个元素是转速测量轴（1轴）上的齿轮齿数， # 第二个元素是与转速测量轴啮合（2轴）的齿轮的齿数， # 当 TransLevel == 2 时，有二对齿轮啮合，3根轴。这个数组有4个元素，第一个元素是转速测量轴（2轴）上的齿轮齿数， # 第二个元素是与转速测量轴啮合（2轴）的齿轮的齿数，第三个元素是2轴与下一级（3轴）啮合的齿轮的齿数， # 第四个元素是3轴与2轴啮合的齿轮的齿数。 #返回： 等角度采样的角度轴坐标，等角度采样的波形，阶次谱的阶次轴坐标，阶次谱的幅值 nRevSpeedSensorShaft = int(len(xTimeSpeedPulse)/nPulsePerRev) # 计算转速测量轴的转速测量信号有多少个整周数 xTimeSpeedPulse = xTimeSpeedPulse[0:0+int(nRevSpeedSensorShaft*nPulsePerRev)] # 整周数以外多余的脉冲信号去掉 xTimeSpeedPerRevSpeedSensorShaft = np.sum(np.reshape(xTimeSpeedPulse,(nRevSpeedSensorShaft,nPulsePerRev)),axis = 1) # 计算转速测量轴每转一周的时间 if TransLevel == 0: xTimeSpeedPerRevCurrentShaft = xTimeSpeedPerRevSpeedSensorShaft xTimeSpeedPerRevCurrentShaft_ave = np.average(xTimeSpeedPerRevCurrentShaft) # 当转速测量轴就是阶次分析的轴时，转速脉冲可以直接使用 fMax_MAX = 1/xTimeSpeedPerRevCurrentShaft_ave*nLine/nRev else: xTemp = np.cumsum(xTimeSpeedPerRevSpeedSensorShaft) # 转速脉冲时间累积 xTemp = np.append(0,xTemp) # 添加起始点 # 当转速测量轴不是阶次分析轴时，传输脉冲信号必须转换到阶次分析轴上 for i in range (TransLevel): up = TransLevelPara[i] # 转速测量轴齿轮齿数 down = TransLevelPara[i+1] # 啮合轴齿轮齿数 dataLength = int(len(xTemp)*up) # 转速测量轴增采样后的数据长度 flipelements = int(len(xTemp)*0.5) # 转速信号加长50%防止长度不够 extention = np.zeros(flipelements) for j in range(flipelements): # 转速信号加长算法，以最后一点为对称中心 extention[j] = xTemp[-1]-(xTemp[-1-j-1]-xTemp[-1]) xTemp1 = np.append(xTemp,extention) # 转速信号加长后的转速脉冲时间累积 xTemp2 = np.zeros(xTemp1.size*(up),xTemp1.dtype) xTemp2[::up] = xTemp1 # 脉冲时间累积信号增采样，增采样率为转速测量轴齿轮的齿数 # plt.figure(30) # plt.plot(xTemp2) #============================================================== delta1 = 0.00000001 delta2 = 0.0000001 freq1 = 0.39/up freq2 = 0.61/up AA = np.matrix([[-4.278e-1,-4.761e-1,0],[-5.941e-1,7.114e-2,0],[-2.66e-3,5.309e-03,0]]) d1 = np.log10(delta1) d2 = np.log10(delta2) D = np.array([1,d1,d1*d1])*AA*np.array([[1],[d2],[d2*d2]]) bb = np.array([11.01217,0.51244]) fK = np.dot(np.array([1.0,(d1-d2)]),bb) df = np.abs(freq2-freq1) N = int(D[0]/df-fK*df+1) if np.remainder(N,2) == 0: N = N+1 upSampleFilterCoefficients = signal.remez(N, [0,0.39,0.61,up/2], [1,0],Hz = up,maxiter = 500) # [freq, response] = signal.freqz(upSampleFilterCoefficients) # ampl = np.abs(response) # plt.figure(31) # plt.plot(freq,ampl) #============================================================== groupDelay = int((len(upSampleFilterCoefficients)+1)/2) xTemp3 = signal.convolve(xTemp2,upSampleFilterCoefficients) xTemp4 = (xTemp3[groupDelay-1:groupDelay-1+dataLength])*up plt.figure(30) plt.plot(xTemp4) xTemp5 = xTemp4[0::down] xTimeSpeedPerRevCurrentShaft = np.diff(xTemp5) xTimeSpeedPerRevCurrentShaft_ave = np.average(xTimeSpeedPerRevCurrentShaft) # 当转速测量轴就是阶次分析的轴时，转速脉冲可以直接使用 fMax_MAX = 1/xTimeSpeedPerRevCurrentShaft_ave*nLine/nRev xTimeSpeedPerRevCurrentShaft_min = np.min(xTimeSpeedPerRevCurrentShaft) # 计算轴转得最快的一周的时间 fMax_MAX = 1/xTimeSpeedPerRevCurrentShaft_min*nLine/nRev # 根据最快的一周的时间得出阶次频谱最高分析频宽的频率 Order = 6 if fMax_MAX <= Fs/2.56: xTimeWave = self.LowPassButter(xTimeWave, Fs, fMax_MAX, Order) # 低通滤波消除大于fMax_MAX的信号，然后才能做阶次分析 xTimeSpeedCurrentShaftCumulative = np.cumsum(xTimeSpeedPerRevCurrentShaft) # 计算轴每转一周所用时间的累积值 nVibDataPoint = nLine*2.56 # 计算时域波形点数 up = nVibDataPoint/nRev # 计算每周等角度采样点数 # down = 1 halfbandWeight = self.halfbandcoefficient() groupDelay = int((len(halfbandWeight)+1)/2) numIteration = int(np.log10(up)/np.log10(2)) temp = np.append(0,xTimeSpeedCurrentShaftCumulative) flipelements = int(len(temp)*0.5) extention = np.zeros(flipelements) for i in range(numIteration): dataLength = int(len(temp)*2) for j in range(flipelements): extention[j] = temp[-1]-(temp[-1-j-1]-temp[-1]) temp1 = np.append(temp,extention) temp2 = np.insert(temp1,np.arange(1,len(temp1),1),0) temp3 = np.append(temp2,0) temp4 = signal.convolve(temp3,halfbandWeight) temp = temp4[groupDelay-1:groupDelay-1+dataLength] xTimeNewSamplePoint = temp flipelements = int(len(xTimeWave)*0.25) extention = np.zeros(flipelements) for j in range(flipelements): extention[j] = xTimeWave[-1]-(xTimeWave[-1-j-1]-xTimeWave[-1]) xTimeWave = np.append(xTimeWave,extention) xTimeSamplePoint = np.arange(0,len(xTimeWave)/Fs,1/Fs) # 计算输入波形的采样时间点 # spl = UnivariateSpline(xTimeSamplePoint,xTimeWave) # Waveform = spl(xTimeNewSamplePoint) fx = interpolate.interp1d(xTimeSamplePoint,xTimeWave,kind= 'cubic') Waveform = fx(xTimeNewSamplePoint) OrderWaveform = Waveform[0:int(nLine*2.56)] ct = 360/(nLine*2.56/nRev) orderTimeLine = np.arange(0,nLine*2.56*ct,ct) fOT = nLine*2.56/nRev fMax = fOT/2.56 fCut = 10*fMax/nLine WinType = 0 AverageType = 0 OverLapType = 0 [f, xSpectrum, xSpectrumPhase] = self.Spectrum(OrderWaveform,nLine,fMax,fCut,WinType, AverageType ,OverLapType) # 计算插值函数 len(OrderWaveform) return [orderTimeLine, OrderWaveform, f, xSpectrum] ############################################################################### ############################################################################### # 内部函数 def halfBandDesign (self, filterLength, transitionBand ): invalidInput = False # check if integer if (np.abs(filterLength - int(filterLength)) > 1e-10): print('halfBandDesign.py: filterLength must be an integer') invalidInput = True # check if too small if (filterLength < 7): print('halfBandDesign.py: filterLength must be larger than 6') invalidInput = True # check if proper length if (np.mod(filterLength+1,4) != 0): print('halfBandDesign.py: filterLength+1 must be divisble by 4') invalidInput = True # check range for transition band if (transitionBand <= 0 or transitionBand >= 0.5): print('halfBandDesign.py: transitionBand must be greater than 0 and less than 0.5') invalidInput = True if (invalidInput): return [] else: # design a half band filter with remez cutoff = 0.25 fPass = cutoff - (transitionBand/2) fStop = cutoff + (transitionBand/2) fVec = [0, fPass, fStop, 0.5] aVec = [1, 0] weights = signal.remez(filterLength,fVec,aVec) # force zero weights zeroWeightIndicesHalf = np.arange(2,(filterLength-1)/2,2,dtype=int) zeroWeightIndicesNegative = np.concatenate((-zeroWeightIndicesHalf[::-1],zeroWeightIndicesHalf)) zeroWeightIndices = zeroWeightIndicesNegative - zeroWeightIndicesNegative[0] + 1 weights[zeroWeightIndices] = 0 return weights ############################################################################### ############################################################################### # def Pk_Search_BW(self): # BW = 0.03 # return BW ############################################################################### ############################################################################### # @staticmethod def gE_Bearing_Defect_HM_Threshhold(self,alert, danger, defectType): gE_Bearing_Defect_HM_alert = np.zeros(5) gE_Bearing_Defect_HM_danger = np.zeros(5) if defectType == 'BPFO' : gE_Bearing_Defect_HM_alert[0] = alert*0.167 gE_Bearing_Defect_HM_alert[1] = alert*0.131 gE_Bearing_Defect_HM_alert[2] = alert*0.098 gE_Bearing_Defect_HM_alert[3] = alert*0.089 gE_Bearing_Defect_HM_alert[4] = alert*0.082 gE_Bearing_Defect_HM_danger[0] = danger*0.167 gE_Bearing_Defect_HM_danger[1] = danger*0.131 gE_Bearing_Defect_HM_danger[2] = danger*0.098 gE_Bearing_Defect_HM_danger[3] = danger*0.089 gE_Bearing_Defect_HM_danger[4] = danger*0.082 return gE_Bearing_Defect_HM_alert, gE_Bearing_Defect_HM_danger if defectType == 'BPFI': gE_Bearing_Defect_HM_alert[0] = alert*0.088 gE_Bearing_Defect_HM_alert[1] = alert*0.064 gE_Bearing_Defect_HM_alert[2] = alert*0.056 gE_Bearing_Defect_HM_alert[3] = alert*0.035 gE_Bearing_Defect_HM_alert[4] = alert*0.024 gE_Bearing_Defect_HM_danger[0] = danger*0.088 gE_Bearing_Defect_HM_danger[1] = danger*0.064 gE_Bearing_Defect_HM_danger[2] = danger*0.056 gE_Bearing_Defect_HM_danger[3] = danger*0.035 gE_Bearing_Defect_HM_danger[4] = danger*0.024 return gE_Bearing_Defect_HM_alert, gE_Bearing_Defect_HM_danger if defectType == 'BSF': gE_Bearing_Defect_HM_alert[0] = alert*0.088 gE_Bearing_Defect_HM_alert[1] = alert*0.064 gE_Bearing_Defect_HM_alert[2] = alert*0.056 gE_Bearing_Defect_HM_alert[3] = alert*0.035 gE_Bearing_Defect_HM_alert[4] = alert*0.024 gE_Bearing_Defect_HM_danger[0] = danger*0.088 gE_Bearing_Defect_HM_danger[1] = danger*0.064 gE_Bearing_Defect_HM_danger[2] = danger*0.056 gE_Bearing_Defect_HM_danger[3] = danger*0.035 gE_Bearing_Defect_HM_danger[4] = danger*0.024 return gE_Bearing_Defect_HM_alert, gE_Bearing_Defect_HM_danger if defectType == 'FTF': gE_Bearing_Defect_HM_alert[0] = alert*0.088 gE_Bearing_Defect_HM_alert[1] = alert*0.064 gE_Bearing_Defect_HM_alert[2] = alert*0.056 gE_Bearing_Defect_HM_alert[3] = alert*0.035 gE_Bearing_Defect_HM_alert[4] = alert*0.024 gE_Bearing_Defect_HM_danger[0] = danger*0.088 gE_Bearing_Defect_HM_danger[1] = danger*0.064 gE_Bearing_Defect_HM_danger[2] = danger*0.056 gE_Bearing_Defect_HM_danger[3] = danger*0.035 gE_Bearing_Defect_HM_danger[4] = danger*0.024 return gE_Bearing_Defect_HM_alert, gE_Bearing_Defect_HM_danger ############################################################################### ############################################################################### # @staticmethod def Vel_Bearing_Defect_HM_Threshhold(self,alert, danger, defectType): Vel_Bearing_Defect_HM_alert = np.zeros(5) Vel_Bearing_Defect_HM_danger = np.zeros(5) if defectType == 'BPFO': Vel_Bearing_Defect_HM_alert[0] = alert*1.0 Vel_Bearing_Defect_HM_alert[1] = alert*1.0 Vel_Bearing_Defect_HM_alert[2] = alert*1.0 Vel_Bearing_Defect_HM_alert[3] = alert*1.0 Vel_Bearing_Defect_HM_danger[0] = danger*1.0 Vel_Bearing_Defect_HM_danger[1] = danger*1.0 Vel_Bearing_Defect_HM_danger[2] = danger*1.0 Vel_Bearing_Defect_HM_danger[3] = danger*1.0 return Vel_Bearing_Defect_HM_alert, Vel_Bearing_Defect_HM_danger if defectType == 'BPFI': Vel_Bearing_Defect_HM_alert[0] = alert*1.0 Vel_Bearing_Defect_HM_alert[1] = alert*1.0 Vel_Bearing_Defect_HM_alert[2] = alert*1.0 Vel_Bearing_Defect_HM_alert[3] = alert*1.0 Vel_Bearing_Defect_HM_danger[0] = danger*1.0 Vel_Bearing_Defect_HM_danger[1] = danger*1.0 Vel_Bearing_Defect_HM_danger[2] = danger*1.0 Vel_Bearing_Defect_HM_danger[3] = danger*1.0 return Vel_Bearing_Defect_HM_alert, Vel_Bearing_Defect_HM_danger if defectType == 'BSF': Vel_Bearing_Defect_HM_alert[0] = alert*1.0 Vel_Bearing_Defect_HM_alert[1] = alert*0.5 Vel_Bearing_Defect_HM_alert[2] = alert*0.5 Vel_Bearing_Defect_HM_alert[3] = alert*0.5 Vel_Bearing_Defect_HM_danger[0] = danger*1.0 Vel_Bearing_Defect_HM_danger[1] = danger*0.5 Vel_Bearing_Defect_HM_danger[2] = danger*0.5 Vel_Bearing_Defect_HM_danger[3] = danger*0.5 return Vel_Bearing_Defect_HM_alert, Vel_Bearing_Defect_HM_danger ############################################################################### ############################################################################### # @staticmethod def Vel_Bearing_RunningSpd_HM_Threshhold(self,alert, danger, defectType): Vel_Bearing_RunningSpd_HM_alert = np.zeros(5) Vel_Bearing_RunningSpd_HM_danger = np.zeros(5) if defectType == 'BSF': Vel_Bearing_RunningSpd_HM_alert[0] = alert*1.0 Vel_Bearing_RunningSpd_HM_alert[1] = alert*0.5 Vel_Bearing_RunningSpd_HM_alert[2] = alert*0.5 Vel_Bearing_RunningSpd_HM_alert[3] = alert*0.5 Vel_Bearing_RunningSpd_HM_alert[4] = alert*0.5 Vel_Bearing_RunningSpd_HM_danger[0] = danger*1.0 Vel_Bearing_RunningSpd_HM_danger[1] = danger*0.5 Vel_Bearing_RunningSpd_HM_danger[2] = danger*0.5 Vel_Bearing_RunningSpd_HM_danger[3] = danger*0.5 Vel_Bearing_RunningSpd_HM_danger[4] = danger*0.5 return Vel_Bearing_RunningSpd_HM_alert, Vel_Bearing_RunningSpd_HM_danger ############################################################################### ############################################################################### # @staticmethod def Vel_Unbalance_HM_Threshold(self,alert, danger): Vel_Unbalance_HM_alert = np.zeros(5) Vel_Unbalance_HM_danger = np.zeros(5) Vel_Unbalance_HM_alert[0] = alert*1.0 Vel_Unbalance_HM_alert[1] = alert*0.5 Vel_Unbalance_HM_alert[2] = alert*0.5 Vel_Unbalance_HM_alert[3] = alert*0.5 Vel_Unbalance_HM_alert[4] = alert*0.5 Vel_Unbalance_HM_danger[0] = danger*1.0 Vel_Unbalance_HM_danger[1] = danger*0.5 Vel_Unbalance_HM_danger[2] = danger*0.5 Vel_Unbalance_HM_danger[3] = danger*0.5 Vel_Unbalance_HM_danger[4] = danger*0.5 return Vel_Unbalance_HM_alert, Vel_Unbalance_HM_danger ############################################################################### ############################################################################### # @staticmethod def Vel_Misalignment_HM_Threshold(self,alert, danger): Vel_Misalignment_HM_alert = np.zeros(5) Vel_Misalignment_HM_danger = np.zeros(5) Vel_Misalignment_HM_alert[0] = alert*1.0 Vel_Misalignment_HM_alert[1] = alert*0.5 Vel_Misalignment_HM_alert[2] = alert*0.5 Vel_Misalignment_HM_alert[3] = alert*0.5 Vel_Misalignment_HM_alert[4] = alert*0.5 Vel_Misalignment_HM_danger[0] = danger*1.0 Vel_Misalignment_HM_danger[1] = danger*0.5 Vel_Misalignment_HM_danger[2] = danger*0.5 Vel_Misalignment_HM_danger[3] = danger*0.5 Vel_Misalignment_HM_danger[4] = danger*0.5 return Vel_Misalignment_HM_alert, Vel_Misalignment_HM_danger ############################################################################### ############################################################################### # @staticmethod def Vel_Looseness_HM_Threshold(self,alert, danger): Vel_Looseness_HM_alert = np.zeros(5) Vel_Looseness_HM_danger = np.zeros(5) Vel_Looseness_HM_alert[0] = alert*1.0 Vel_Looseness_HM_alert[1] = alert*0.5 Vel_Looseness_HM_alert[2] = alert*0.5 Vel_Looseness_HM_alert[3] = alert*0.5 Vel_Looseness_HM_alert[4] = alert*0.5 Vel_Looseness_HM_danger[0] = danger*1.0 Vel_Looseness_HM_danger[1] = danger*0.5 Vel_Looseness_HM_danger[2] = danger*0.5 Vel_Looseness_HM_danger[3] = danger*0.5 Vel_Looseness_HM_danger[4] = danger*0.5 return Vel_Looseness_HM_alert, Vel_Looseness_HM_danger ############################################################################### ############################################################################### # @staticmethod def gE_2XLF_HM_Threshhold(self,alert, danger): gE_2XLF_HM_alert = np.zeros(5) gE_2XLF_HM_danger = np.zeros(5) gE_2XLF_HM_alert[0] = alert*0.167 gE_2XLF_HM_alert[1] = alert*0.131 gE_2XLF_HM_alert[2] = alert*0.098 gE_2XLF_HM_alert[3] = alert*0.089 gE_2XLF_HM_alert[4] = alert*0.082 gE_2XLF_HM_danger[0] = danger*0.167 gE_2XLF_HM_danger[1] = danger*0.131 gE_2XLF_HM_danger[2] = danger*0.098 gE_2XLF_HM_danger[3] = danger*0.089 gE_2XLF_HM_danger[4] = danger*0.082 return gE_2XLF_HM_alert, gE_2XLF_HM_danger ############################################################################### ############################################################################### # @staticmethod def Vel_2XLF_HM_Threshhold(self,alert, danger): Vel_2XLF_HM_alert = np.zeros(1) Vel_2XLF_HM_danger = np.zeros(1) Vel_2XLF_HM_alert[0] = alert*1.0 Vel_2XLF_HM_danger[0] = danger*1.0 return Vel_2XLF_HM_alert, Vel_2XLF_HM_danger ############################################################################### ############################################################################### # @staticmethod def Acc_GMF_HM_Threshhold(self,alert, danger): Acc_GMF_HM_alert = np.zeros(3) Acc_GMF_HM_danger = np.zeros(3) Acc_GMF_HM_alert[0] = alert*1.0 Acc_GMF_HM_alert[1] = alert*1.0 Acc_GMF_HM_alert[2] = alert*1.0 Acc_GMF_HM_danger[0] = danger*1.0 Acc_GMF_HM_danger[1] = danger*1.0 Acc_GMF_HM_danger[2] = danger*1.0 return Acc_GMF_HM_alert, Acc_GMF_HM_danger ############################################################################### ############################################################################### # @staticmethod def Vel_GMF_HM_Threshhold(self,alert, danger): Vel_GMF_HM_alert = np.zeros(3) Vel_GMF_HM_danger = np.zeros(3) Vel_GMF_HM_alert[0] = alert*1.0 Vel_GMF_HM_alert[1] = alert*1.0 Vel_GMF_HM_alert[2] = alert*1.0 Vel_GMF_HM_danger[0] = danger*1.0 Vel_GMF_HM_danger[1] = danger*1.0 Vel_GMF_HM_danger[2] = danger*1.0 return Vel_GMF_HM_alert, Vel_GMF_HM_danger ############################################################################### ############################################################################### # @staticmethod def gE_GMF_HM_Threshhold(self,alert, danger): gE_GMF_HM_alert = np.zeros(3) gE_GMF_HM_danger = np.zeros(3) gE_GMF_HM_alert[0] = alert*0.167 gE_GMF_HM_alert[1] = alert*0.131 gE_GMF_HM_alert[2] = alert*0.098 gE_GMF_HM_danger[0] = danger*0.167 gE_GMF_HM_danger[1] = danger*0.131 gE_GMF_HM_danger[2] = danger*0.098 return gE_GMF_HM_alert, gE_GMF_HM_danger ############################################################################### ############################################################################### # @staticmethod def gE_GTIF_HM_Threshhold(self,alert, danger): gE_GTIF_HM_alert = np.zeros(5) gE_GTIF_HM_danger = np.zeros(5) gE_GTIF_HM_alert[0] = alert*0.167 gE_GTIF_HM_alert[1] = alert*0.131 gE_GTIF_HM_alert[2] = alert*0.098 gE_GTIF_HM_alert[3] = alert*0.089 gE_GTIF_HM_alert[4] = alert*0.082 gE_GTIF_HM_danger[0] = danger*0.167 gE_GTIF_HM_danger[1] = danger*0.131 gE_GTIF_HM_danger[2] = danger*0.098 gE_GTIF_HM_danger[3] = danger*0.089 gE_GTIF_HM_danger[4] = danger*0.082 return gE_GTIF_HM_alert, gE_GTIF_HM_danger ############################################################################### ############################################################################### # @staticmethod def Vel_VPF_HM_Threshhold(self,alert, danger): Vel_VPF_HM_alert = np.zeros(2) Vel_VPF_HM_danger = np.zeros(2) Vel_VPF_HM_alert[0] = alert*1.0 Vel_VPF_HM_alert[1] = alert*1.0 Vel_VPF_HM_danger[0] = danger*1.0 Vel_VPF_HM_danger[1] = danger*1.0 return Vel_VPF_HM_alert, Vel_VPF_HM_danger ############################################################################### ############################################################################### # @staticmethod def Vel_BTF_HM_Threshhold(self,alert, danger): Vel_BTF_HM_alert = np.zeros(2) Vel_BTF_HM_danger = np.zeros(2) Vel_BTF_HM_alert[0] = alert*1.0 Vel_BTF_HM_alert[1] = alert*1.0 Vel_BTF_HM_danger[0] = danger*1.0 Vel_BTF_HM_danger[1] = danger*1.0 return Vel_BTF_HM_alert, Vel_BTF_HM_danger ############################################################################### ############################################################################### # @staticmethod def Vel_SHF_HM_Threshhold(self,alert, danger): Vel_SHF_HM_alert = np.zeros(1) Vel_SHF_HM_danger = np.zeros(1) Vel_SHF_HM_alert[0] = alert*1.0 Vel_SHF_HM_danger[0] = danger*1.0 return Vel_SHF_HM_alert, Vel_SHF_HM_danger ############################################################################### ############################################################################### # @staticmethod def Vel_RBPF_HM_Threshhold(self,alert, danger): Vel_RBPF_HM_alert = np.zeros(1) Vel_RBPF_HM_danger = np.zeros(1) Vel_RBPF_HM_alert[0] = alert*1.0 Vel_RBPF_HM_danger[0] = danger*1.0 return Vel_RBPF_HM_alert, Vel_RBPF_HM_danger ############################################################################### ############################################################################### # @staticmethod def Acc_RBPF_HM_Threshhold(self,alert, danger): Acc_RBPF_HM_alert = np.zeros(1) Acc_RBPF_HM_danger = np.zeros(1) Acc_RBPF_HM_alert[0] = alert*1.0 Acc_RBPF_HM_danger[0] = danger*1.0 return Acc_RBPF_HM_alert, Acc_RBPF_HM_danger ############################################################################### ############################################################################### # @staticmethod def Bearing_Fault_Diag(self,gE_Bearing_Defect_HM, Vel_Bearing_Defect_HM, Vel_RunningSpd_HM, gE_alert, gE_danger, Vel_alert, Vel_danger, defectType): BrDefDetection = False BrDefSev = 0 gE_Bearing_Defect_HM_alert = np.zeros(5) gE_Bearing_Defect_HM_danger = np.zeros(5) Vel_Bearing_Defect_HM_alert = np.zeros(5) Vel_Bearing_Defect_HM_danger = np.zeros(5) Vel_Bearing_RunningSpd_HM_alert = np.zeros(5) Vel_Bearing_RunningSpd_HM_danger = np.zeros(5) ############################################################################### # Check bearing outer race and inner race defect if (defectType == 'BPFO') or (defectType == 'BPFI'): [gE_Bearing_Defect_HM_alert, gE_Bearing_Defect_HM_danger] = self.gE_Bearing_Defect_HM_Threshhold(gE_alert, gE_danger, defectType) # [Vel_Bearing_RunningSpd_HM_alert, Vel_Bearing_RunningSpd_HM_danger] = self.Vel_Bearing_RunningSpd_HM_Threshhold(Vel_alert, Vel_danger, defectType) [Vel_Bearing_Defect_HM_alert, Vel_Bearing_Defect_HM_danger] = self.Vel_Bearing_Defect_HM_Threshhold(Vel_alert, Vel_danger, defectType) # Bearing outer race frequency, severe, confirmed by velocity spectrum, CF = 8, Severity = 10 # AND # OR # ge_BPFO_1 >= Alert # ge_BPFO_2 >= Alert # OR # Vel_BPFO_1 >= Alert # Vel_BPFO_2 >= Alert # # above logics apply to BPFI as well if (gE_Bearing_Defect_HM[0] >= gE_Bearing_Defect_HM_alert[0] or gE_Bearing_Defect_HM[1] >= gE_Bearing_Defect_HM_alert[1]) and \ (Vel_Bearing_Defect_HM[0] >= Vel_Bearing_Defect_HM_alert[0] or Vel_Bearing_Defect_HM[1] >= Vel_Bearing_Defect_HM_alert[1]): BrDefDetection = True BrDefSev = 8 return BrDefDetection, BrDefSev # Bearing outer race frequency, severe, no confirming by velocity spectrum, CF = 6, Severity = 10 # AND # OR # AND # ge_BPFO_1 >= Danger # ge_BPFO_2 >= Alert # AND # gE_BPFO_1 >= Danger # gE_BPFO_3 >= Alert # AND # gE_BPFO_1 >= Danger # gE_BPFO_4 >= Alert # NOR # Vel_BPFO_1 >= Alert # Vel_BPFO_2 >= Alert # # above logics apply to BPFI as well else: if ((gE_Bearing_Defect_HM[0] >= gE_Bearing_Defect_HM_danger[0] and gE_Bearing_Defect_HM[1] >= gE_Bearing_Defect_HM_alert[1]) or \ (gE_Bearing_Defect_HM[0] >= gE_Bearing_Defect_HM_danger[0] and gE_Bearing_Defect_HM[2] >= gE_Bearing_Defect_HM_alert[2]) or \ (gE_Bearing_Defect_HM[0] >= gE_Bearing_Defect_HM_danger[0] and gE_Bearing_Defect_HM[3] >= gE_Bearing_Defect_HM_alert[3])) and \ (not(Vel_Bearing_Defect_HM[0] >= Vel_Bearing_Defect_HM_alert[0] or Vel_Bearing_Defect_HM[1] >= Vel_Bearing_Defect_HM_alert[1])): BrDefDetection = True BrDefSev = 6 return BrDefDetection, BrDefSev # Bearing outer race frequency, severe, not confirmed by gE spectrum, CF = 6, Severity = 10 # OR # Vel_BPFO_1 >= Alert # Vel_BPFO_2 >= Alert # # above logics apply to BPFI as well else: if (Vel_Bearing_Defect_HM[0] >= Vel_Bearing_Defect_HM_alert[0] or Vel_Bearing_Defect_HM[1] >= Vel_Bearing_Defect_HM_alert[1]): BrDefDetection = True BrDefSev = 6 return BrDefDetection, BrDefSev # Bearing outer race frequency, moderate, several defect peaks, no confirming by velocity spectrum, CF = 5, Severity = 5 # AND # OR # AND # ge_BPFO_1 >= Alert # ge_BPFO_2 >= Alert # AND # gE_BPFO_1 >= Alert # gE_BPFO_3 >= Alert # AND # gE_BPFO_1 >= Alert # gE_BPFO_4 >= Alert # NOR # Vel_BPFO_1 >= Alert # Vel_BPFO_2 >= Alert # # above logics apply to BPFI as well else: if ((gE_Bearing_Defect_HM[0] >= gE_Bearing_Defect_HM_alert[0] and gE_Bearing_Defect_HM[1] >= gE_Bearing_Defect_HM_alert[1]) or \ (gE_Bearing_Defect_HM[0] >= gE_Bearing_Defect_HM_alert[0] and gE_Bearing_Defect_HM[2] >= gE_Bearing_Defect_HM_alert[2]) or \ (gE_Bearing_Defect_HM[0] >= gE_Bearing_Defect_HM_alert[0] and gE_Bearing_Defect_HM[3] >= gE_Bearing_Defect_HM_alert[3])) and \ (not(Vel_Bearing_Defect_HM[0] >= Vel_Bearing_Defect_HM_alert[0] or Vel_Bearing_Defect_HM[1] >= Vel_Bearing_Defect_HM_alert[1])): BrDefDetection = True BrDefSev = 3 return BrDefDetection, BrDefSev # Bearing outer race frequency, suspect, only 1x defect peak evaluated, CF = 3, Severity = 5 # ge_BPFO_1 >= Alert # # above logics apply to BPFI as well else: if (gE_Bearing_Defect_HM[0] >= gE_Bearing_Defect_HM_alert[0]): BrDefDetection = True BrDefSev = 1 return BrDefDetection, BrDefSev return BrDefDetection, BrDefSev ############################################################################### # Check bearing rolling element defect if defectType == 'BSF': [gE_Bearing_Defect_HM_alert, gE_Bearing_Defect_HM_danger] = self.gE_Bearing_Defect_HM_Threshhold(gE_alert, gE_danger, defectType) [Vel_Bearing_RunningSpd_HM_alert, Vel_Bearing_RunningSpd_HM_danger] = self.Vel_Bearing_RunningSpd_HM_Threshhold(Vel_alert, Vel_danger, defectType) [Vel_Bearing_Defect_HM_alert, Vel_Bearing_Defect_HM_danger] = self.Vel_Bearing_Defect_HM_Threshhold(Vel_alert, Vel_danger, defectType) # print(Vel_Bearing_RunningSpd_HM_alert, Vel_Bearing_RunningSpd_HM_danger) # Bearing rolling element frequency, severe, CF = 9, Severity = 10 # AND # AND # Vel_BSF_2 >= Alert # Vel_BSF_4 >= Alert # NOR # Vel_RS_2 >= Alert # Vel_RS_3 >= Alert # Vel_RS_4 >= Alert # Vel_RS_5 >= Alert if (Vel_Bearing_Defect_HM[1] >= Vel_Bearing_Defect_HM_alert[1] and Vel_Bearing_Defect_HM[3] >= Vel_Bearing_Defect_HM_alert[3]) and \ (not(Vel_RunningSpd_HM[1] >= Vel_Bearing_RunningSpd_HM_alert[1] or Vel_RunningSpd_HM[2] >= Vel_Bearing_RunningSpd_HM_alert[2] or \ Vel_RunningSpd_HM[3] >= Vel_Bearing_RunningSpd_HM_alert[3] or Vel_RunningSpd_HM[4] >= Vel_Bearing_RunningSpd_HM_alert[4])): BrDefDetection = True BrDefSev = 9 return BrDefDetection, BrDefSev # Bearing rolling element frequency, moderate, CF = 9, Severity = 5 # AND # AND # gE_BSF_1 >= Alert # NOR # Vel_BSF_2 >= Alert # Vel_BSF_4 >= Alert # Vel_RS_2 >= Alert # Vel_RS_3 >= Alert # Vel_RS_4 >= Alert # Vel_RS_5 >= Alert else: if gE_Bearing_Defect_HM[0] >= gE_Bearing_Defect_HM_alert[0] and \ (not(Vel_Bearing_Defect_HM[1] >= Vel_Bearing_Defect_HM_alert[1] or Vel_Bearing_Defect_HM[3] >= Vel_Bearing_Defect_HM_alert[3] or \ Vel_RunningSpd_HM[1] >= Vel_Bearing_RunningSpd_HM_alert[1] or Vel_RunningSpd_HM[2] >= Vel_Bearing_RunningSpd_HM_alert[2] or \ Vel_RunningSpd_HM[3] >= Vel_Bearing_RunningSpd_HM_alert[3] or Vel_RunningSpd_HM[4] >= Vel_Bearing_RunningSpd_HM_alert[4])): BrDefDetection = True BrDefSev = 5 return BrDefDetection, BrDefSev # Bearing rolling element frequency, suspect, CF = 3, Severity = 5 # AND # OR # gE_BSF_1 = Alert # gE_BSF_2 = Alert # gE_BSF_4 = Alert # NOR # gE_BSF_1 = Danger # gE_BSF_2 = Danger # gE_BSF_4 = Danger # Vel_BSF_2 >= Alert # Vel_BSF_4 >= Alert # Vel_RS_2 >= Alert # Vel_RS_3 >= Alert # Vel_RS_4 >= Alert # Vel_RS_5 >= Alert else: if ((gE_Bearing_Defect_HM[0] >= gE_Bearing_Defect_HM_alert[0] and gE_Bearing_Defect_HM[0] < gE_Bearing_Defect_HM_danger[0]) or \ (gE_Bearing_Defect_HM[1] >= gE_Bearing_Defect_HM_alert[1] and gE_Bearing_Defect_HM[1] < gE_Bearing_Defect_HM_danger[1]) or \ (gE_Bearing_Defect_HM[3] >= gE_Bearing_Defect_HM_alert[3] and gE_Bearing_Defect_HM[3] < gE_Bearing_Defect_HM_danger[3] )) and \ (not(gE_Bearing_Defect_HM[0] >= gE_Bearing_Defect_HM_danger[0] or gE_Bearing_Defect_HM[1] >= gE_Bearing_Defect_HM_danger[1] or gE_Bearing_Defect_HM[3] >= gE_Bearing_Defect_HM_danger[3] or \ Vel_Bearing_Defect_HM[0] >= Vel_Bearing_Defect_HM_alert[0] or Vel_Bearing_Defect_HM[1] >= Vel_Bearing_Defect_HM_alert[1] or \ Vel_RunningSpd_HM[1] >= Vel_Bearing_RunningSpd_HM_alert[1] or Vel_RunningSpd_HM[2] >= Vel_Bearing_RunningSpd_HM_alert[2] or \ Vel_RunningSpd_HM[3] >= Vel_Bearing_RunningSpd_HM_alert[3] or Vel_RunningSpd_HM[4] >= Vel_Bearing_RunningSpd_HM_alert[4])): BrDefDetection = True BrDefSev = 1 return BrDefDetection, BrDefSev return BrDefDetection, BrDefSev ############################################################################### # Check bearing cage defect if defectType == 'FTF': [gE_Bearing_Defect_HM_alert, gE_Bearing_Defect_HM_danger] = self.gE_Bearing_Defect_HM_Threshhold(gE_alert, gE_danger, defectType) # Bearing cage frequency, severe, CF = 6, Severity = 10 # OR # gE_CG_1 >= Danger # gE_CG_2 >= Danger if gE_Bearing_Defect_HM[0] >= gE_Bearing_Defect_HM_danger[0] or gE_Bearing_Defect_HM[1] >= gE_Bearing_Defect_HM_danger[1]: BrDefDetection = True BrDefSev = 6 return BrDefDetection, BrDefSev # Bearing cage frequency, moderate, CF = 5, Severity = 5 # AND # OR # gE_CG_1 >= Alert # gE_CG_2 >= Alert # NOR # gE_CG_1 = Danger # gE_CG_2 = Danger else: if (gE_Bearing_Defect_HM[0] >= gE_Bearing_Defect_HM_alert[0] or gE_Bearing_Defect_HM[1] >= gE_Bearing_Defect_HM_alert[1]) and \ (not(gE_Bearing_Defect_HM[0] >= gE_Bearing_Defect_HM_danger[0] or gE_Bearing_Defect_HM[1] >= gE_Bearing_Defect_HM_danger[1])): BrDefDetection = True BrDefSev = 3 return BrDefDetection, BrDefSev return BrDefDetection, BrDefSev ############################################################################### ############################################################################### # @staticmethod def RunningSpd_Symptom_Rule(self,Vel_RunningSpd_HM, Vel_RunningSpd_HM_alert, Vel_RunningSpd_HM_danger): RPM1_MOD = False RPM1_SEV = False RPM2_MOD = False RPM2_SEV = False RPM345_MOD = False RPM345_SEV = False # 1x RPM, Moderate, CF=8, Severity = 5 # AND # Vel_RS_1 = Alert if Vel_RunningSpd_HM[0] >= Vel_RunningSpd_HM_alert[0] and Vel_RunningSpd_HM[0] < Vel_RunningSpd_HM_danger[0]: RPM1_MOD = True # 1x RPM, Severe, CF=8, Severity = 10 # AND # Vel_RS_1 = Danger if Vel_RunningSpd_HM[0] >= Vel_RunningSpd_HM_danger[0]: RPM1_SEV = True # 2x RPM, Moderate, CF=8, Severity = 5 # AND # Vel_RS_2 = Alert if Vel_RunningSpd_HM[1] >= Vel_RunningSpd_HM_alert[1] and Vel_RunningSpd_HM[1] < Vel_RunningSpd_HM_danger[1]: RPM2_MOD = True # 2x RPM, Severe, CF=8, Severity = 10 # AND # Vel_RS_2 = Danger if Vel_RunningSpd_HM[1] >= Vel_RunningSpd_HM_danger[1]: RPM2_SEV = True # 3x, 4x, 5x RPM, moderate, CF=8, Severity = 5 # AND # OR # Vel_RS_3 = Alert # Vel_RS_4 = Alert # Vel_RS_5 = Alert # NOR # Vel_RS_3 = Danger # Vel_RS_4 = Danger # Vel_RS_5 = Danger if (Vel_RunningSpd_HM[2] >= Vel_RunningSpd_HM_alert[2] and Vel_RunningSpd_HM[2] < Vel_RunningSpd_HM_danger[2] or \ Vel_RunningSpd_HM[3] >= Vel_RunningSpd_HM_alert[3] and Vel_RunningSpd_HM[3] < Vel_RunningSpd_HM_danger[3] or \ Vel_RunningSpd_HM[4] >= Vel_RunningSpd_HM_alert[4] and Vel_RunningSpd_HM[4] < Vel_RunningSpd_HM_danger[4]) and \ (not(Vel_RunningSpd_HM[2] >= Vel_RunningSpd_HM_danger[2] or Vel_RunningSpd_HM[3] >= Vel_RunningSpd_HM_danger[3] or \ Vel_RunningSpd_HM[4] >= Vel_RunningSpd_HM_danger[4])): RPM345_MOD = True # 3x, 4x, 5x RPM, Sever, CF=8, Severity = 10 # OR # Vel_RS_3 = Danger # Vel_RS_4 = Danger # Vel_RS_5 = Danger if Vel_RunningSpd_HM[2] >= Vel_RunningSpd_HM_danger[2] or Vel_RunningSpd_HM[3] >= Vel_RunningSpd_HM_danger[3] or \ Vel_RunningSpd_HM[4] >= Vel_RunningSpd_HM_danger[4]: RPM345_SEV = True return RPM1_MOD, RPM1_SEV, RPM2_MOD, RPM2_SEV, RPM345_MOD, RPM345_SEV ############################################################################### ############################################################################### # @staticmethod # Note: The alert and danger threshold values (in mm/s RMS) are set according to ISO vibration standard for flexible mounting in group 3 def RunningSpd_Fault_Diag(self,Vel_RunningSpd_HM, Vel_alert = 7.1, Vel_danger = 11, defectType = 'Unbalance'): RPM1_MOD = False RPM1_SEV = False RPM2_MOD = False RPM2_SEV = False RPM345_MOD = False RPM345_SEV = False UnbalanceDetection = False UnbalanceSev = 0 MisalignmentDetction = False MisalignmentSev = 0 LoosenessDetection = False LoosenessSev = 0 Vel_Unbalance_HM_alert = np.zeros(5) Vel_Unbalance_HM_danger = np.zeros(5) Vel_Misalignment_HM_alert = np.zeros(5) Vel_Misalignment_HM_danger = np.zeros(5) Vel_Looseness_HM_alert = np.zeros(5) Vel_Looseness_HM_danger = np.zeros(5) ############################################################################### # Check Unbalance if defectType == 'Unbalance': [Vel_Unbalance_HM_alert, Vel_Unbalance_HM_danger] = self.Vel_Unbalance_HM_Threshold(Vel_alert, Vel_danger) [RPM1_MOD, RPM1_SEV, RPM2_MOD, RPM2_SEV, RPM345_MOD, RPM345_SEV] = self.RunningSpd_Symptom_Rule(Vel_RunningSpd_HM, Vel_Unbalance_HM_alert, Vel_Unbalance_HM_danger) # print('RPM1_MOD,RPM1_SEV, RPM2_MOD, RPM2_SEV, RPM345_MOD, RPM345_SEV',' ',RPM1_MOD,' ', RPM1_SEV,' ', RPM2_MOD, ' ',RPM2_SEV, ' ',RPM345_MOD, ' ',RPM345_SEV) # Unbalance, severe, rather confident, 1x RPM peak dominates, CF=8, Severity = 10 # AND # RPM1_SEV == True # NOR # RPM2_MOD == True # RPM2_SEV == True # RPM345_MOD == True # RPM345_SEV == True if RPM1_SEV == True and (not (RPM2_MOD == True or RPM2_SEV == True or RPM345_MOD == True or RPM345_SEV == True )): UnbalanceDetection = True UnbalanceSev = 8 else: # Unbalance, severe, somewhat confident 1x RPM peak dominates but low level of 2xRPM is possible, CF=6, Severity = 10 # AND # RPM1_SEV == True # NOR # RPM2_SEV == True # RPM345_MOD == True # RPM345_SEV == True if RPM1_SEV == True and (not(RPM2_SEV == True or RPM345_MOD == True or RPM345_SEV == True )): UnbalanceDetection = True UnbalanceSev = 6 else: # Unbalance, moderate, rather confident since 1x RPM peak dominates, CF=8, Severity = 5 # AND # RPM1_MOD == True # NOR # RPM2_MOD == True # RPM2_SEV == True # RPM345_MOD == True # RPM345_SEV == True if RPM1_MOD == True and (not(RPM2_MOD == True or RPM2_SEV == True or RPM345_MOD == True or RPM345_SEV == True )): UnbalanceDetection = True UnbalanceSev = 5 return UnbalanceDetection, UnbalanceSev ############################################################################### # Check Misalignment if defectType == 'Misalignment': [Vel_Misalignment_HM_alert, Vel_Misalignment_HM_danger] = self.Vel_Misalignment_HM_Threshold(Vel_alert, Vel_danger) [RPM1_MOD, RPM1_SEV, RPM2_MOD, RPM2_SEV, RPM345_MOD, RPM345_SEV] = self.RunningSpd_Symptom_Rule(Vel_RunningSpd_HM, Vel_Misalignment_HM_alert, Vel_Misalignment_HM_danger) # Misalignment, moderate CF = 8, Severity = 10 # AND # RPM2_SEV == True # OR # RPM1_MOD == True # RPM1_SEV == True # RPM345_MOD == True # NOR # RPM345_SEV == True if RPM2_SEV == True and (RPM1_MOD == True or RPM1_SEV == True or RPM345_MOD == True) and (not(RPM1_SEV == True)): MisalignmentDetction = True MisalignmentSev = 8 # Misalignment, moderate CF = 6, Severity 10 # AND # RPM2_SEV == True # NOR # RPM1_MOD == True # RPM1_SEV == True # RPM345_MOD == True # RPM345_SEV == True else: if RPM2_SEV == True and (not(RPM1_MOD == True or RPM1_SEV == True or RPM345_MOD == True or RPM345_SEV == True)): MisalignmentDetction = True MisalignmentSev = 6 # Misalignment, moderate CF = 8, Severity = 5 # AND # RPM2_MOD == True # OR # RPM1_MOD == True # RPM345_MOD == True # NOR # RPM1_SEV == True # RPM345_SEV == True else: if RPM2_MOD == True and (RPM1_MOD == True or RPM345_MOD == True) and (not(RPM1_SEV == True or RPM345_SEV == True)): MisalignmentDetction = True MisalignmentSev = 5 # Misalignment, moderate CF = 6, Severity = 5 # AND # RPM2_MOD == True # NOR # RPM1_MOD == True # RPM1_SEV == True # RPM345_MOD == True # RPM345_SEV == True else: if RPM2_MOD == True and (not (RPM1_MOD == True or RPM1_SEV == True or RPM345_MOD == True or RPM345_SEV == True)): MisalignmentDetction = True MisalignmentSev = 4 return MisalignmentDetction, MisalignmentSev ############################################################################### # Check Looseness if defectType == 'Looseness': [Vel_Looseness_HM_alert, Vel_Looseness_HM_danger] = self.Vel_Looseness_HM_Threshold(Vel_alert, Vel_danger) [RPM1_MOD, RPM1_SEV, RPM2_MOD, RPM2_SEV, RPM345_MOD, RPM345_SEV] = self.RunningSpd_Symptom_Rule(Vel_RunningSpd_HM, Vel_Looseness_HM_alert, Vel_Looseness_HM_danger) # Mechanical looseness, sever CF = 8, Severity = 10 # OR # AND # RPM1_SEV == True # OR # RPM345_MOD == True # RPM345_SEV == True # AND # RPM2_SEV == True # OR # RPM345_MOD == True # RPM345_SEV == True # AND # RPM345_SEV == True # OR # RPM1_MOD == True # RPM1_SEV == True # AND # RPM345_SEV == True # OR # RPM2_MOD == True # RPM2_SEV == True if (RPM1_SEV == True and (RPM345_MOD == True or RPM345_SEV == True )) or (RPM2_SEV == True and (RPM345_MOD == True or RPM345_SEV == True )) or \ (RPM345_SEV == True and (RPM1_MOD == True or RPM1_SEV == True )) or (RPM345_SEV == True and (RPM2_MOD == True or RPM2_SEV == True )): LoosenessDetection = True LoosenessSev = 8 # Mechanical looseness, sever CF = 6, Severity = 10 # AND # RPM345_SEV == True # NOR # RPM1_MOD == True # RPM1_SEV == True # RPM2_MOD == True # RPM2_SEV == True else: if RPM345_SEV == True and (not( RPM1_MOD == True or RPM1_SEV == True or RPM2_MOD == True or RPM2_SEV == True )): LoosenessDetection = True LoosenessSev = 6 # Mechanical looseness, moderate CF = 8, Severity = 5 # AND # OR # AND # RPM1_MOD == True # RPM345_MOD == True # AND # RPM2_MOD == True # RPM345_MOD == True # NOR # RPM1_SEV == True # RPM2_SEV == True # RPM345_SEV == True else: if ((RPM1_MOD == True and RPM345_MOD == True )or(RPM2_MOD == True and RPM345_MOD == True))and\ (not(RPM1_SEV == True or RPM2_SEV == True or RPM345_SEV == True)): LoosenessDetection = True LoosenessSev = 5 # Mechanical looseness, moderate CF = 6, Severity = 5 # AND # RPM345_MOD == True # NOR # RPM1_MOD == True # RPM1_SEV == True # RPM2_MOD == True # RPM2_SEV == True else: if RPM345_MOD == True and (not( RPM1_MOD == True or RPM1_SEV == True or RPM2_MOD == True or \ RPM2_SEV == True )): LoosenessDetection = True LoosenessSev = 4 return LoosenessDetection, LoosenessSev ############################################################################### ############################################################################### # @staticmethod def GMF_GTIF_Symptom_Rule(self,Vel_GMF_HM, Acc_GMF_HM, gE_GMF_HM, gE_GTIF_HM, Vel_GMF_HM_alert, Vel_GMF_HM_danger, \ Acc_GMF_HM_alert, Acc_GMF_HM_danger, gE_GMF_HM_alert, gE_GMF_HM_danger, gE_GTIF_HM_alert, \ gE_GTIF_HM_danger): GMF_MOD = False GMF_SEV = False GTIF = False # Gear mesh frequency and harmonics, moderate, CF=8, Severity = 5 # AND # OR # Vel_GMF_1 = Alert # Vel_GMF_2 = Alert # Vel_GMF_3 = Alert # Acc_GMF_1 = Alert # Acc_GMF_2 = Alert # Acc_GMF_3 = Alert # gE_GMF_1 = Alert # gE_GMF_2 = Alert # gE_GMF_3 = Alert # NOR # Vel_GMF_1 = Danger # Vel_GMF_2 = Danger # Vel_GMF_3 = Danger # Acc_GMF_1 = Danger # Acc_GMF_2 = Danger # Acc_GMF_3 = Danger # gE_GMF_1 = Danger # gE_GMF_2 = Danger # gE_GMF_3 = Danger if (Vel_GMF_HM[0] >= Vel_GMF_HM_alert[0] and Vel_GMF_HM[0] < Vel_GMF_HM_danger[0] or \ Vel_GMF_HM[1] >= Vel_GMF_HM_alert[1] and Vel_GMF_HM[1] < Vel_GMF_HM_danger[1] or \ Vel_GMF_HM[2] >= Vel_GMF_HM_alert[2] and Vel_GMF_HM[2] < Vel_GMF_HM_danger[2] or \ Acc_GMF_HM[0] >= Acc_GMF_HM_alert[0] and Acc_GMF_HM[0] < Acc_GMF_HM_danger[0] or \ Acc_GMF_HM[1] >= Acc_GMF_HM_alert[1] and Acc_GMF_HM[1] < Acc_GMF_HM_danger[1] or \ Acc_GMF_HM[2] >= Acc_GMF_HM_alert[2] and Acc_GMF_HM[2] < Acc_GMF_HM_danger[2] or \ gE_GMF_HM[0] >= gE_GMF_HM_alert[0] and gE_GMF_HM[0] < gE_GMF_HM_danger[0] or \ gE_GMF_HM[1] >= gE_GMF_HM_alert[1] and gE_GMF_HM[1] < gE_GMF_HM_danger[1] or \ gE_GMF_HM[2] >= gE_GMF_HM_alert[2] and gE_GMF_HM[2] < gE_GMF_HM_danger[2]) \ and \ (not(Vel_GMF_HM[0] >= Vel_GMF_HM_danger[0] or \ Vel_GMF_HM[1] >= Vel_GMF_HM_danger[1] or \ Vel_GMF_HM[2] >= Vel_GMF_HM_danger[2] or \ Acc_GMF_HM[0] >= Acc_GMF_HM_danger[0] or \ Acc_GMF_HM[1] >= Acc_GMF_HM_danger[1] or \ Acc_GMF_HM[2] >= Acc_GMF_HM_danger[2] or \ gE_GMF_HM[0] >= gE_GMF_HM_danger[0] or \ gE_GMF_HM[1] >= gE_GMF_HM_danger[1] or \ gE_GMF_HM[2] >= gE_GMF_HM_danger[2])): GMF_MOD = True # Gear mesh frequency and harmonics, sever, CF=8, Severity = 10 # OR # Vel_GMF_1 = Danger # Vel_GMF_2 = Danger # Vel_GMF_3 = Danger # Acc_GMF_1 = Danger # Acc_GMF_2 = Danger # Acc_GMF_3 = Danger # gE_GMF_1 = Danger # gE_GMF_2 = Danger # gE_GMF_3 = Danger if (Vel_GMF_HM[0] >= Vel_GMF_HM_danger[0] or \ Vel_GMF_HM[1] >= Vel_GMF_HM_danger[1] or \ Vel_GMF_HM[2] >= Vel_GMF_HM_danger[2] or \ Acc_GMF_HM[0] >= Acc_GMF_HM_danger[0] or \ Acc_GMF_HM[1] >= Acc_GMF_HM_danger[1] or \ Acc_GMF_HM[2] >= Acc_GMF_HM_danger[2] or \ gE_GMF_HM[0] >= gE_GMF_HM_danger[0] or \ gE_GMF_HM[1] >= gE_GMF_HM_danger[1] or \ gE_GMF_HM[2] >= gE_GMF_HM_danger[2]): GMF_SEV = True # Gear tooth impact frequency and harmonics, CF=8, Severity = 10 # OR # gE_GTIF_1 >= Alert # gE_GTIF_2 >= Alert # gE_GTIF_3 >= Alert # gE_GTIF_4 >= Alert # gE_GTIF_5 >= Alert if (gE_GTIF_HM[0] >= gE_GTIF_HM_alert[0] or \ gE_GTIF_HM[1] >= gE_GTIF_HM_alert[1] or \ gE_GTIF_HM[2] >= gE_GTIF_HM_alert[2] or \ gE_GTIF_HM[3] >= gE_GTIF_HM_alert[3] or \ gE_GTIF_HM[4] >= gE_GTIF_HM_alert[4]): GTIF = True return GMF_MOD, GMF_SEV, GTIF ############################################################################### ############################################################################### # @staticmethod def Gear_Fault_Diag(self,Acc_GMF_HM, Vel_GMF_HM, gE_GMF_HM, gE_GTIF_HM, Acc_GMF_alert_gPk = 0.299, Acc_GMF_danger_gPk = 0.344, \ Vel_GMF_alert_mmsRMS = 4.0, Vel_GMF_danger_mmsRMS = 8.0, gE_GMF_alert_gE = 3.0, \ gE_GMF_danger_gE = 6.0, gE_GTIF_alert_gE = 3.0, gE_GTIF_danger_gE = 6.0): GMF_MOD = False GMF_SEV = False GTIF = False Acc_GMF_HM_alert =[] Acc_GMF_HM_danger = [] Vel_GMF_HM_alert =[] Vel_GMF_HM_danger = [] gE_GMF_HM_alert = [] gE_GMF_HM_danger = [] gE_GTIF_HM_alert = [] gE_GTIF_HM_danger = [] GearToothWearDetection = False GearToothWearSev = 0 GearToothCrackBrokenDetection = False GearToothCrackBrokenSev = 0 [Acc_GMF_HM_alert, Acc_GMF_HM_danger] = self.Acc_GMF_HM_Threshhold(Acc_GMF_alert_gPk, Acc_GMF_danger_gPk) [Vel_GMF_HM_alert, Vel_GMF_HM_danger] = self.Vel_GMF_HM_Threshhold(Vel_GMF_alert_mmsRMS, Vel_GMF_danger_mmsRMS) [gE_GMF_HM_alert, gE_GMF_HM_danger] = self.gE_GMF_HM_Threshhold(gE_GMF_alert_gE, gE_GMF_danger_gE) [gE_GTIF_HM_alert, gE_GTIF_HM_danger] = self.gE_GTIF_HM_Threshhold(gE_GTIF_alert_gE, gE_GTIF_danger_gE) [GMF_MOD, GMF_SEV, GTIF] = self.GMF_GTIF_Symptom_Rule(Vel_GMF_HM, Acc_GMF_HM, gE_GMF_HM, gE_GTIF_HM, Vel_GMF_HM_alert, \ Vel_GMF_HM_danger, Acc_GMF_HM_alert, Acc_GMF_HM_danger, \ gE_GMF_HM_alert, gE_GMF_HM_danger, gE_GTIF_HM_alert, gE_GTIF_HM_danger) # Gear wear, misalignment, or eccentricity, sever, CF = 8, severity = 10 # AND # GMF_SEV = True if GMF_SEV == True: GearToothWearDetection = True GearToothWearSev = 10 # Gear wear, misalignment, or eccentricity, moderate, CF = 8, severity = 5 # AND # GMF_MOD = True else: if GMF_MOD == True: GearToothWearDetection = True GearToothWearSev = 5 # Gear tooth broken or cracked, CF = 8, severity = 5 # AND # GTIF = True # NOR # GMF_MOD = True # GMF_SEV = True if GTIF == True and (not(GMF_MOD == True or GMF_SEV == True)): GearToothCrackBrokenDetection = True GearToothCrackBrokenSev = 5 return GearToothWearDetection, GearToothWearSev, GearToothCrackBrokenDetection, GearToothCrackBrokenSev ############################################################################### ############################################################################### # @staticmethod # Note: the default alert and danger threshold values (in mm/s RMS for velocity and g peak for acceleration) are set according a DSS AC motor model, it can be adjusted if necessary def RBPF_Symptom_Rule(self,Vel_RBPF_HM, Acc_RBPF_HM, Vel_RBPF_HM_alert, Vel_RBPF_HM_danger, Acc_RBPF_HM_alert, \ Acc_RBPF_HM_danger): RBPF_MOD = False RBPF_SEV = False # Rotorbar pass frequency, sever, CF=8, severity = 10 # OR # Vel_RBPF_1 == Danger # Acc_RBPF_1 == Danger if (Vel_RBPF_HM[0] >= Vel_RBPF_HM_danger[0] or Acc_RBPF_HM[0] >= Acc_RBPF_HM_danger[0]): RBPF_SEV = True # Rotorbar pass frequency, moderate, CF=8, severity = 5 # OR # Vel_RBPF_1 == Alert # Acc_RBPF_1 == Alert else: if (Vel_RBPF_HM[0] >= Vel_RBPF_HM_alert[0] or Acc_RBPF_HM[0] >= Acc_RBPF_HM_alert[0]) and Vel_RBPF_HM[0] < Vel_RBPF_HM_danger[0] and Acc_RBPF_HM[0] < Acc_RBPF_HM_danger[0]: RBPF_MOD = True return RBPF_MOD, RBPF_SEV ############################################################################### ############################################################################### # @staticmethod # Note: the default alert and danger threshold values (in mm/s RMS for velocity and g peak for acceleration) are set according a DSS AC motor model, it can be adjusted if necessary def Cracked_Rotorbar_Diag(self,Vel_RBPF_HM, Acc_RBPF_HM, Vel_RBPF_alert_mmsRMS = 0.064*25.4, \ Vel_RBPF_danger_mmsRMS = 0.127*25.4, Acc_RBPF_alert_gPk = 2.5, Acc_RBPF_danger_gPk = 3.75): RBPF_MOD = False RBPF_SEV = False Vel_RBPF_HM_alert = [] Vel_RBPF_HM_danger = [] Acc_RBPF_HM_alert = [] Acc_RBPF_HM_danger = [] CrackedRotorbarDetection = False CrackedRotorbarSev = 0 [Acc_RBPF_HM_alert, Acc_RBPF_HM_danger] = self.Acc_RBPF_HM_Threshhold(Acc_RBPF_alert_gPk, Acc_RBPF_danger_gPk) [Vel_RBPF_HM_alert,Vel_RBPF_HM_danger] = self.Vel_RBPF_HM_Threshhold(Vel_RBPF_alert_mmsRMS, Vel_RBPF_danger_mmsRMS) [RBPF_MOD, RBPF_SEV] = self.RBPF_Symptom_Rule(Vel_RBPF_HM, Acc_RBPF_HM, Vel_RBPF_HM_alert, Vel_RBPF_HM_danger, Acc_RBPF_HM_alert, Acc_RBPF_HM_danger) # Cracked rotorbar, sever CF = 8, severity = 10 # AND # RBPF_SEV = True if RBPF_SEV == True: CrackedRotorbarDetection = True CrackedRotorbarSev = 8 # Cracked rotorbar, moderate CF = 8, severity = 5 # AND # RBPF_MOD = True # NOR # RBPF_SEV = True if RBPF_MOD == True and (not(RBPF_SEV == True)): CrackedRotorbarDetection = True CrackedRotorbarSev = 5 return CrackedRotorbarDetection, CrackedRotorbarSev ############################################################################### ############################################################################### # @staticmethod # Note: the default velocity alert and danger threshold (in mm/s RMS for velocity) is set according a DSS AC motor model, it can be adjusted if necessary def LF2X_Symptom_Rule(self,Vel_2XLF_HM, gE_2XLF_HM, Vel_2XLF_HM_alert, Vel_2XLF_HM_danger, gE_2XLF_HM_alert, gE_2XLF_HM_danger): LF2X_MOD = False LF2X_SEV = False # 2X line frequency, moderate, CF=8, severity = 5 # AND # Vel_2XLF_1 = alert if Vel_2XLF_HM[0] >= Vel_2XLF_HM_alert[0] and Vel_2XLF_HM[0] < Vel_2XLF_HM_danger[0]: LF2X_MOD = True # 2X line frequency, moderate, several defect peaks, no confirming velocity, CF=8, severity = 5 # AND # OR # AND # gE_2XLF_1 = alert # gE_2XLF_2 = alert # AND # gE_2XLF_1 = alert # gE_2XLF_3 = alert # AND # gE_2XLF_1 = alert # gE_2XLF_4 = alert # AND # gE_2XLF_1 = alert # gE_2XLF_5 = alert # NOR # Vel_2XLF_1 >= alert else: if (((gE_2XLF_HM[0] >= gE_2XLF_HM_alert[0] and gE_2XLF_HM[0] < gE_2XLF_HM_danger[0]) and (gE_2XLF_HM[1] >= gE_2XLF_HM_alert[1] and gE_2XLF_HM[1] < gE_2XLF_HM_danger[1])) or \ ((gE_2XLF_HM[0] >= gE_2XLF_HM_alert[0] and gE_2XLF_HM[0] < gE_2XLF_HM_danger[0]) and (gE_2XLF_HM[2] >= gE_2XLF_HM_alert[2] and gE_2XLF_HM[2] < gE_2XLF_HM_danger[2])) or \ ((gE_2XLF_HM[0] >= gE_2XLF_HM_alert[0] and gE_2XLF_HM[0] < gE_2XLF_HM_danger[0]) and (gE_2XLF_HM[3] >= gE_2XLF_HM_alert[3] and gE_2XLF_HM[3] < gE_2XLF_HM_danger[3])) or \ ((gE_2XLF_HM[0] >= gE_2XLF_HM_alert[0] and gE_2XLF_HM[0] < gE_2XLF_HM_danger[0]) and (gE_2XLF_HM[4] >= gE_2XLF_HM_alert[4] and gE_2XLF_HM[4] < gE_2XLF_HM_danger[4]))) and \ (not(Vel_2XLF_HM[0] >= Vel_2XLF_HM_alert[0])): LF2X_MOD = True # 2X line frequency, severe, CF=8, severity = 10 # AND # Vel_2XLF_1 = Danger if Vel_2XLF_HM[0] >= Vel_2XLF_HM_danger[0]: LF2X_SEV = True # 2X line frequency, severe, several defect peaks, no confirming velocity, CF=8, severity = 10 # AND # OR # AND # gE_2XLF_1 >= danger # gE_2XLF_2 >= alert # AND # gE_2XLF_1 >= danger # gE_2XLF_3 >= alert # AND # gE_2XLF_1 >= danger # gE_2XLF_4 >= alert # AND # gE_2XLF_1 >= danger # gE_2XLF_5 >= alert # NOR # Vel_2XLF_1 >= alert else: if (((gE_2XLF_HM[0] >= gE_2XLF_HM_danger[0]) and (gE_2XLF_HM[1] >= gE_2XLF_HM_alert[1])) or \ ((gE_2XLF_HM[0] >= gE_2XLF_HM_danger[0]) and (gE_2XLF_HM[2] >= gE_2XLF_HM_alert[2])) or \ ((gE_2XLF_HM[0] >= gE_2XLF_HM_danger[0]) and (gE_2XLF_HM[3] >= gE_2XLF_HM_alert[3])) or \ ((gE_2XLF_HM[0] >= gE_2XLF_HM_danger[0]) and (gE_2XLF_HM[4] >= gE_2XLF_HM_alert[4]))) and \ (not(Vel_2XLF_HM[0] < Vel_2XLF_HM_alert[0])): LF2X_SEV = True # 2X line frequency, severe, confirmed by velocity spectrum, CF=8, severity = 10 # AND # OR # gE_2XLF_1 >= Alert # gE_2XLF_2 >= Alert # OR # Vel_2XLF_1 >= Alert else: if (gE_2XLF_HM[0] >= gE_2XLF_HM_alert[0] or gE_2XLF_HM[1] >= gE_2XLF_HM_alert[1]) and \ Vel_2XLF_HM[0] >= Vel_2XLF_HM_alert[0]: LF2X_SEV = True return LF2X_MOD, LF2X_SEV ############################################################################### ############################################################################### # @staticmethod def Variable_AirGap_Diag(self,Vel_2XLF_HM, gE_2XLF_HM, Vel_2XLF_alert = 0.064*25.4, Vel_2XLF_danger = 0.127*25.4, gE_2XLF_alert = 3, gE_2XLF_danger = 10): LF2X_MOD = False LF2X_SEV = False gE_2XLF_HM_alert = [] gE_2XLF_HM_danger = [] Vel_2XLF_HM_alert = [] Vel_2XLF_HM_danger = [] VariableAirgapDetection = False VariableAirgapSev = 0 [gE_2XLF_HM_alert,gE_2XLF_HM_danger] = self.gE_2XLF_HM_Threshhold(gE_2XLF_alert, gE_2XLF_danger) [Vel_2XLF_HM_alert,Vel_2XLF_HM_danger] = self.Vel_2XLF_HM_Threshhold(Vel_2XLF_alert, Vel_2XLF_danger) [LF2X_MOD, LF2X_SEV] = self.LF2X_Symptom_Rule(Vel_2XLF_HM, gE_2XLF_HM, Vel_2XLF_HM_alert, Vel_2XLF_HM_danger,gE_2XLF_HM_alert, gE_2XLF_HM_danger) # Variable air gap, sever CF = 8, severity = 10 # AND # LF2X_SEV = True if LF2X_SEV == True: VariableAirgapDetection = True VariableAirgapSev = 8 # Variable air gap, moderate CF = 8, severity = 5 # AND # LF2X_MOD = True # NOR # LF2X_SEV = True else: if LF2X_MOD == True and LF2X_SEV == False: VariableAirgapDetection = True VariableAirgapSev = 5 return VariableAirgapDetection, VariableAirgapSev ############################################################################### ############################################################################### def VPF_Symptom_Rule(self,Vel_VPF_HM, Vel_VPF_HM_alert, Vel_VPF_HM_danger): VAN_VPFH_MOD = False VAN_VPFH_SEV = False # Vane pass frequency and harmonics (VPFH), sever, CF=8, Severity = 10 # OR # Vel_VPF_1 = Danger # Vel_VPF_2 = Danger if Vel_VPF_HM[0] >= Vel_VPF_HM_danger[0] or Vel_VPF_HM[1] >= Vel_VPF_HM_danger[1] : VAN_VPFH_SEV = True # Vane pass frequency and harmonics (VPFH), moderate, CF=8, Severity = 5 # AND # OR # Vel_VPF_1 = Alert # Vel_VPF_2 = Alert # NOR # Vel_VPF_1 = Danger # Vel_VPF_2 = Danger if (Vel_VPF_HM[0] >= Vel_VPF_HM_alert[0] and Vel_VPF_HM[0] < Vel_VPF_HM_danger[0] or\ Vel_VPF_HM[1] >= Vel_VPF_HM_alert[1] and Vel_VPF_HM[1] < Vel_VPF_HM_danger[1]) and\ (not(Vel_VPF_HM[0] >= Vel_VPF_HM_danger[0] or Vel_VPF_HM[1] >= Vel_VPF_HM_danger[1])): VAN_VPFH_MOD = True return VAN_VPFH_MOD, VAN_VPFH_SEV ############################################################################### ############################################################################### # @staticmethod def Vane_WearDamage_Diag(self,Vel_VPF_HM, Vel_VPF_alert_mmsRMS = 0.159*25.4, \ Vel_VPF_danger_mmsRMS = 0.318*25.4): VAN_VPFH_MOD = False VAN_VPFH_SEV = False Vel_VPF_HM_alert = [] Vel_VPF_HM_danger = [] VaneWearDamageDetection = False VaneWearDamageSev = 0 [Vel_VPF_HM_alert,Vel_VPF_HM_danger] = self.Vel_VPF_HM_Threshhold(Vel_VPF_alert_mmsRMS, Vel_VPF_danger_mmsRMS) [VAN_VPFH_MOD, VAN_VPFH_SEV] = self.VPF_Symptom_Rule(Vel_VPF_HM, Vel_VPF_HM_alert, Vel_VPF_HM_danger) # Vane wear or damage, severe CF = 8, severity = 10 # OR # VAN_VPFH_SEV= True if VAN_VPFH_SEV == True: VaneWearDamageDetection = True VaneWearDamageSev = 10 # Vane wear or damage, moderate CF = 6, severity = 5 # OR # VAN_VPFH_MOD= True else: if VAN_VPFH_MOD == True: VaneWearDamageDetection = True VaneWearDamageSev = 5 return VaneWearDamageDetection, VaneWearDamageSev ############################################################################### ############################################################################### # 巴特沃斯带通滤波 # xTimeWave是时域波形数组， fs是采样频率，单位为Hz # lowcut 为低截止频率，单位为Hz # highcut 为高截止频率，单位为Hz # Order为滤波器阶数，通常为偶数，2，4，6，8一般不超过8 def BandPassButter(self,xTimeWave, lowcut, highcut, fs, order): low = lowcut * 2 / fs high = highcut * 2 / fs sos = signal.butter(order, [low, high], 'bandpass', output='sos') filtered = signal.sosfilt(sos, xTimeWave) return filtered