Browse Source

更新热失控预警,使用M-K趋势检验算法

qingfeng 3 years ago
parent
commit
888d163918

+ 5 - 5
LIB/MIDDLE/CellStateEstimation/BatSafetyWarning/V1_0_1/CBMSBatInterShort.py

@@ -120,10 +120,10 @@ class BatInterShort():
     def _celldeltsoc_get(self,cellvolt_list,dict_baltime,capacity): 
         cellsoc=[]
         celldeltsoc=[]
-        for j in range(self.param.CellVoltNums):    #获取每个电芯电压对应的SOC值
+        for j in range(self.param.CellVoltNums):   #获取每个电芯电压对应的SOC值
             cellvolt=cellvolt_list[j]
             ocv_soc=np.interp(cellvolt,self.param.LookTab_OCV,self.param.LookTab_SOC)
-            if j in dict_baltime.keys():
+            if j+1 in dict_baltime.keys():
                 ocv_soc=ocv_soc+dict_baltime[j+1]*self.param.BalCurrent/(capacity*3600)   #补偿均衡电流
             else:
                 pass
@@ -337,7 +337,7 @@ class BatInterShort():
 
         #容量初始化
         if self.df_soh.empty:
-            batsoh=self.df_bms.loc[0,'SOH[%]']
+            batsoh=(self.df_bms.loc[0,'SOH[%]'])
             capacity=self.param.Capacity*batsoh/100
         else:
             batsoh=self.df_soh.loc[len(self.df_soh)-1,'soh']
@@ -540,7 +540,7 @@ class BatInterShort():
 
         #容量初始化
         if self.df_soh.empty:
-            batsoh=self.df_bms.loc[0,'SOH[%]']
+            batsoh=(self.df_bms.loc[0,'SOH[%]'])
             capacity=self.param.Capacity*batsoh/100
         else:
             batsoh=self.df_soh.loc[len(self.df_soh)-1,'soh']
@@ -639,7 +639,7 @@ class BatInterShort():
                         pass                
                 elif standingtime>3600*12:
                     standingtime=0
-                    cellvolt_now=np.array(self._avgvolt_get(i))
+                    cellvolt_now=self._avgvolt_get(i)
                     if not cellvolt_now.empty:
                         cellvolt_min=min(cellvolt_now)
                         cellvolt_max=max(cellvolt_now)

+ 92 - 19
LIB/MIDDLE/CellStateEstimation/BatSafetyWarning/V1_0_1/CBMSSafetyWarning.py

@@ -3,6 +3,7 @@ import numpy as np
 import datetime
 import time
 from matplotlib import pyplot as plt
+import pymannkendall as mk
 from LIB.MIDDLE.CellStateEstimation.Common.V1_0_1 import BatParam
 
 class SafetyWarning:
@@ -31,6 +32,7 @@ class SafetyWarning:
         df_res=pd.DataFrame(columns=['start_time', 'end_time', 'product_id', 'code', 'level', 'info','advice'])
 
         time_now=datetime.datetime.now()
+        time_now=time_now.strftime('%Y-%m-%d %H:%M:%S')
         time_sp='0000-00-00 00:00:00'
 
         #参数初始化.......................................
@@ -41,6 +43,11 @@ class SafetyWarning:
         R2_list=[]
         voltsigmafault_list=[]
         uniformfault_list=[]
+        mk_trend_list=[]
+        mk_p_list=[]
+        mk_z_list=[]
+        mk_Tau_list=[]
+        mk_slope_list=[]
 
         if not self.df_short.empty:
             short_current=self.df_short['short_current']
@@ -103,34 +110,63 @@ class SafetyWarning:
             #电压变化率及电压离群度.................................................................................
             if not self.OutLineVol_Rate.empty and VoltChange.iloc[-1]['time']*36000>18*3600 and len(VoltChange)>5:
 
-                volt3sigma=np.array(Volt_3Sigma[volt_column[i]].map(lambda x:eval(x)))
-                volt3sigma_sum=np.sum(volt3sigma<-3)
-                #电压变化率
                 VoltChange[volt_column[i]]=VoltChange[volt_column[i]].map(lambda x:eval(x))
                 y=VoltChange[volt_column[i]]
-                a1,b1=np.polyfit(VoltChange['time'].tolist(),y.tolist(),1)
-                y1=a1*VoltChange['time']+b1
-                y_mean=y.mean()
-                R2=1-(np.sum((y1-y)**2))/(np.sum((y-y_mean)**2))
-                R2_list.append(R2)
+                volt3sigma=np.array(Volt_3Sigma[volt_column[i]].map(lambda x:eval(x)))
+                volt3sigma_sum=np.sum(volt3sigma<-3)
+                # #电压变化率
+                # a1,b1=np.polyfit(VoltChange['time'].tolist(),y.tolist(),1)
+                # y1=a1*VoltChange['time']+b1
+                # y_mean=y.mean()
+                # R2=1-(np.sum((y1-y)**2))/(np.sum((y-y_mean)**2))
+                # R2_list.append(R2)
                 
-                volt_rate.append(a1)
-                # plt.plot(xtime1,y1,label=a1)
+                # volt_rate.append(a1)
+                # plt.plot(xtime1,y1,label='单体电压'+str(i+1))
+                # plt.xlabel('时间', fontsize=25)
+                # plt.ylabel('SOC差', fontsize=25)
+                # plt.xticks(fontsize=20)
+                # plt.yticks(fontsize=20)
                 # plt.scatter( xtime1,VoltChange[volt_column[i]],marker='o')
-                # plt.legend(loc='best')
-                # plt.title('单体电压'+str(i+1))
+                # plt.legend(bbox_to_anchor=(1, 0), loc=3, fontsize=16)
+                # plt.title(self.celltype)
                 # plt.rcParams['font.sans-serif']=['SimHei'] #用来正常显示中文标签
                 # plt.rcParams['axes.unicode_minus']=False #用来正常显示负号  
-                # plt.show()
+                # # plt.show()
                 
                 if volt3sigma_sum>len(volt3sigma)/2:
                     voltsigmafault=1
                 else:
                     voltsigmafault=0
+                voltsigmafault_list.append(voltsigmafault)
+                
+                #mana-kendell趋势检验
+                mk_res=mk.regional_test(np.array(y)) 
+                mk_trend_list.append(mk_res.trend)
+                mk_p_list.append(mk_res.p)
+                mk_z_list.append(mk_res.z)
+                mk_Tau_list.append(mk_res.Tau)
+                mk_slope_list.append(mk_res.slope)
+                """
+                trend:趋势;
+                h:有无趋势;
+                p:趋势的显著水平,越小趋势越明显;
+                z:检验统计量,正代表随时间增大趋势,负代表随时间减小趋势;
+                Tau:反映两个序列的相关性,接近1的值表示强烈的正相关,接近-1的值表示强烈的负相关;
+                s:如果S是一个正数,那么后一部分的观测值相比之前的观测值会趋向于变大;如果S是一个负数,那么后一部分的观测值相比之前的观测值会趋向于变小
+                slope:趋势斜率
+                """
+                # print('单体电压{}:\n'.format(i+1), mk_res)
+
             else:
                 volt_rate.append(0)
                 R2_list.append(0)
-            voltsigmafault_list.append(voltsigmafault)
+                voltsigmafault_list.append(0)
+                mk_trend_list.append(0)
+                mk_p_list.append(0)
+                mk_z_list.append(0)
+                mk_Tau_list.append(0)
+                mk_slope_list.append(0)
 
             #电芯SOC排名判断.............................................................................
             if not self.df_uniform.empty:
@@ -151,13 +187,24 @@ class SafetyWarning:
                 faultinfo='电芯{}发生热失控安全预警'.format(i+1)
                 faultadvice='1联系用户远离电池,立刻召回电池'
                 df_res.loc[0]=[time_now, time_sp, self.sn, faultcode, faultlv, faultinfo, faultadvice]
+                break
             else:
                 pass       
         
-        #电池电压变化率热失控预警...............................................................................
+        #电池电压变化率离群度计算...............................................................................
         volt_rate_std=np.std(volt_rate)
         volt_rate_mean=np.mean(volt_rate)
         volt_rate_3sigma=(np.array(volt_rate)-volt_rate_mean)/volt_rate_std
+        
+        #mk离群度计算
+        mk_slope_std=np.std(mk_slope_list)
+        mk_slope_mean=np.mean(mk_slope_list)
+        mk_slope_3sigma=(np.array(mk_slope_list)-mk_slope_mean)/mk_slope_std   
+
+        mk_z_std=np.std(mk_z_list)
+        mk_z_mean=np.mean(mk_z_list)
+        mk_z_3sigma=(np.array(mk_z_list)-mk_z_mean)/mk_z_std  
+
         if not self.df_soh.empty and self.celltype<50:
             cellsoh=eval(self.df_soh.loc[0,'cellsoh'])
             cellsoh_std=np.std(cellsoh)
@@ -165,8 +212,35 @@ class SafetyWarning:
             cellsoh_3sigma=((np.array(cellsoh)-cellsoh_mean)/cellsoh_std)
         else:
             cellsoh_3sigma=[0]*self.param.CellVoltNums
-        for i in range(len(volt_rate)):
-            if volt_rate[i]<self.param.TrwVoltRate and volt_rate_3sigma[i]<-3 and abs(cellsoh_3sigma[i])<2.5 and  voltsigmafault_list[i]==1:
+        
+        #电压/SOC变化率
+        # for i in range(len(volt_rate)):
+        #     if volt_rate[i]<self.param.TrwVoltRate and volt_rate_3sigma[i]<-3 and abs(cellsoh_3sigma[i])<2.5 and  voltsigmafault_list[i]==1:
+        #         faultcode=110
+        #         faultlv=4
+        #         faultinfo='电芯{}发生热失控安全预警'.format(i+1)
+        #         faultadvice='2联系用户远离电池,立刻召回电池'
+        #         df_res.loc[0]=[time_now, time_sp, self.sn, faultcode, faultlv, faultinfo, faultadvice]
+        #     else:
+        #         pass
+
+        #mana-kendall趋势检测
+        for i in range(len(mk_p_list)):
+            #适用动态工况判断
+            if mk_trend_list[i]=='decreasing' and mk_p_list[i]<self.param.mk_p and mk_z_list[i]<self.param.mk_z and mk_Tau_list[i]<self.param.mk_Tau and mk_slope_3sigma[i]<-3 and mk_slope_list[i]<self.param.mk_slope and abs(cellsoh_3sigma[i])<2.5:
+                faultcode=110
+                faultlv=4
+                faultinfo='电芯{}发生热失控安全预警'.format(i+1)
+                faultadvice='2联系用户远离电池,立刻召回电池'
+                df_res.loc[0]=[time_now, time_sp, self.sn, faultcode, faultlv, faultinfo, faultadvice]
+            #适用静态工况判断
+            elif self.celltype<=50 and mk_trend_list[i]=='decreasing' and mk_p_list[i]<self.param.mk_p and mk_z_list[i]<-6 and (mk_Tau_list[i]<-0.9 or mk_z_3sigma[i]<-3) and mk_slope_3sigma[i]<-3.5 and mk_slope_list[i]<-0.025:
+                faultcode=110
+                faultlv=4
+                faultinfo='电芯{}发生热失控安全预警'.format(i+1)
+                faultadvice='2联系用户远离电池,立刻召回电池'
+                df_res.loc[0]=[time_now, time_sp, self.sn, faultcode, faultlv, faultinfo, faultadvice]
+            elif self.celltype>50 and mk_trend_list[i]=='decreasing' and mk_p_list[i]<self.param.mk_p and mk_z_list[i]<-6 and (mk_Tau_list[i]<-0.9 or mk_z_3sigma[i]<-3) and mk_slope_3sigma[i]<-3.5 and mk_slope_list[i]<-0.25:
                 faultcode=110
                 faultlv=4
                 faultinfo='电芯{}发生热失控安全预警'.format(i+1)
@@ -174,7 +248,6 @@ class SafetyWarning:
                 df_res.loc[0]=[time_now, time_sp, self.sn, faultcode, faultlv, faultinfo, faultadvice]
             else:
                 pass
-        
-        # plt.show()
 
+        # plt.show()
         return df_res

+ 3 - 4
LIB/MIDDLE/CellStateEstimation/BatSafetyWarning/V1_0_1/VoltStray.py

@@ -79,8 +79,8 @@ def main(sn,df_bms,celltype):
     df_bms['time']=pd.to_datetime(df_bms['时间戳'], format='%Y-%m-%d %H:%M:%S')
     volt_column = ['单体电压'+str(i) for i in range(1,param.CellVoltNums+1)]
     columns=['time']+volt_column
-    df_bms=df_bms[df_bms['SOC[%]']>10]
-    # df_standing=df_bms[(df_bms['PackCrnt']>-0.1) & (df_bms['PackCrnt']<0.1)]
+    df_bms=df_bms[(df_bms['SOC[%]']>10)]
+    df_bms=df_bms[(df_bms['PackCrnt']<1)]
     # df_chrg=df_bms[(df_bms['PackCrnt']<-1)]
 
     #电压/SOC变化率计算
@@ -99,7 +99,6 @@ def main(sn,df_bms,celltype):
 
         VolChng = cal_volt_change(df_soc,volt_column)
     else:
-        df_bms=df_bms[df_bms['SOC[%]']<95]
         # df_bms=df_bms[(df_bms['PackCrnt']>-0.1) & (df_bms['PackCrnt']<0.1)]
         df_ori = df_bms[columns]
         df = df_ori.drop_duplicates(subset=['time']) # 删除时间相同的数据
@@ -117,7 +116,7 @@ def main(sn,df_bms,celltype):
 
     OutLineVol=DataFrame(columns=['time','sn','VolOl_Uni','VolChng_Uni'])
 
-    #静置电压变化率和离群度计算
+    #电压变化率和离群度计算
     if len(VolChng)>5 and len(VolSigma)>5:
         VolChng['time'] = time_list1
         VolChng= VolChng.set_index('time')

+ 16 - 8
LIB/MIDDLE/CellStateEstimation/Common/V1_0_1/BatParam.py

@@ -29,6 +29,12 @@ class BatParam:
         self.SohLow=70
         self.SohDiff=15
 
+        #mana-kendall趋势检验参数
+        self.mk_p=0.001
+        self.mk_z=-5
+        self.mk_Tau=-0.7
+        self.mk_slope=-0.1
+
         self.OcvWeight_Temp=[-30,-20,-10,0,10,20,30,40,50]
         self.OcvWeight_StandingTime=[0,500,600,1200,1800,3600,7200,10800]
         self.OcvWeight            =[[0,0,  0,  0,    0,   0.1,0.3, 1],
@@ -91,7 +97,7 @@ class BatParam:
             self.CellTempDiffLv1=10
             self.CellTempDiffLv2=15
 
-            self.TrwVoltRate=-1
+            # self.TrwVoltRate=-1
             
             self.DifVolGap = 3
             self.CellOVlmt=5
@@ -136,7 +142,7 @@ class BatParam:
             self.CellTempDiffLv1=10
             self.CellTempDiffLv2=15
 
-            self.TrwVoltRate=-1
+            # self.TrwVoltRate=-1
             
             self.DifVolGap = 3
             self.CellOVlmt=5
@@ -181,7 +187,7 @@ class BatParam:
             self.CellTempDiffLv1=10
             self.CellTempDiffLv2=15  
 
-            self.TrwVoltRate=-1
+            # self.TrwVoltRate=-1
             
             self.DifVolGap = 3
             self.CellOVlmt=5
@@ -201,8 +207,8 @@ class BatParam:
             self.PackCrntDec=-1
             self.BalCurrent=0.015
 
-            self.LookTab_SOC = [-5,   0,	    5,	    10,	    15,	    20,	    25,	    30,	    35,	    40,	    45,	    50,	    55,	    60,	    65,	    70,	    75,	    80,	    85,	    90,	    95,	    100,   105,     110]
-            self.LookTab_OCV = [3.0,  3.152, 	3.397, 	3.438, 	3.481, 	3.523, 	3.560, 	3.586, 	3.604, 	3.620, 	3.638, 	3.661, 	3.693, 	3.748, 	3.803, 	3.853, 	3.903, 	3.953, 	4.006, 	4.063, 	4.121, 	4.183, 4.253,  4.50]
+            self.LookTab_SOC = [0,	    5,	    10,	    15,	    20,	    25,	    30,	    35,	    40,	    45,	    50,	    55,	    60,	    65,	    70,	    75,	    80,	    85,	    90,	    95,	    100,   105,     110]
+            self.LookTab_OCV = [3.152, 	3.397, 	3.438, 	3.481, 	3.523, 	3.560, 	3.586, 	3.604, 	3.620, 	3.638, 	3.661, 	3.693, 	3.748, 	3.803, 	3.853, 	3.903, 	3.953, 	4.006, 	4.063, 	4.121, 	4.183, 4.253,  4.50]
 
 
             self.CellOvLv1=4.3
@@ -229,7 +235,7 @@ class BatParam:
             self.CellTempDiffLv1=10
             self.CellTempDiffLv2=15   
 
-            self.TrwVoltRate=-1    
+            # self.TrwVoltRate=-1    
             
             self.DifVolGap = 3
             self.CellOVlmt=5
@@ -285,7 +291,8 @@ class BatParam:
             self.CellTempDiffLv1=10
             self.CellTempDiffLv2=15   
 
-            self.TrwVoltRate=-8
+            # self.TrwVoltRate=-8
+            self.mk_slope=-0.3
             
             
             self.DifVolGap = 3
@@ -339,7 +346,8 @@ class BatParam:
             self.CellTempDiffLv1=28
             self.CellTempDiffLv2=32   
 
-            self.TrwVoltRate=-8  
+            # self.TrwVoltRate=-8 
+            self.mk_slope=-0.3
             
                                    
             self.DifVolGap = 3