HealthAssessor.py 17 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412
  1. from functools import lru_cache
  2. from typing import Dict, List
  3. import numpy as np
  4. import pandas as pd
  5. from sklearn.neighbors import BallTree
  6. from app.logger import logger
  7. class HealthAssessor:
  8. def __init__(self):
  9. self.subsystem_config = {
  10. # 发电机
  11. 'generator': {
  12. # 双馈
  13. 'dfig': {
  14. 'fixed': ['generator_winding1_temperature', 'generator_winding2_temperature',
  15. 'generator_winding3_temperature', 'generatordrive_end_bearing_temperature',
  16. 'generatornon_drive_end_bearing_temperature'],
  17. },
  18. # 直驱
  19. 'direct': {
  20. 'fixed': ['generator_winding1_temperature', 'generator_winding2_temperature',
  21. 'generator_winding3_temperature', 'main_bearing_temperature'],
  22. }
  23. },
  24. # 机舱系统
  25. 'nacelle': {
  26. 'fixed': ['front_back_vibration_of_the_cabin', 'side_to_side_vibration_of_the_cabin',
  27. 'cabin_position', 'cabin_temperature'],
  28. },
  29. # 电网环境
  30. 'grid': {
  31. 'fixed': ['reactive_power', 'active_power', 'grid_a_phase_current',
  32. 'grid_b_phase_current', 'grid_c_phase_current'],
  33. },
  34. # 传动系统
  35. 'drive_train': {
  36. 'fixed': ['main_bearing_temperature'],
  37. 'keywords': [
  38. {'include': ['gearbox', 'temperature'], 'exclude': [], 'min_count': 2},
  39. ]
  40. }
  41. }
  42. # 嵌入源代码的MSET实现
  43. self.mset = self._create_mset_core()
  44. def _create_mset_core(self):
  45. """创建MSET核心计算模块"""
  46. class MSETCore:
  47. def __init__(self):
  48. self.matrixD = None
  49. self.normalDataBallTree = None
  50. self.healthyResidual = None
  51. def calcSimilarity(self, x, y):
  52. """优化后的相似度计算"""
  53. diff = np.array(x) - np.array(y)
  54. return 1 / (1 + np.sqrt(np.sum(diff ** 2)))
  55. def genDLMatrix(self, trainDataset, dataSize4D=60, dataSize4L=5):
  56. """优化矩阵生成过程"""
  57. m, n = trainDataset.shape
  58. # 快速选择极值点
  59. min_indices = np.argmin(trainDataset, axis=0)
  60. max_indices = np.argmax(trainDataset, axis=0)
  61. unique_indices = np.unique(np.concatenate([min_indices, max_indices]))
  62. self.matrixD = trainDataset[unique_indices].copy()
  63. # 快速填充剩余点
  64. remaining_indices = np.setdiff1d(np.arange(m), unique_indices)
  65. np.random.shuffle(remaining_indices)
  66. needed = max(0, dataSize4D - len(unique_indices))
  67. if needed > 0:
  68. self.matrixD = np.vstack([self.matrixD, trainDataset[remaining_indices[:needed]]])
  69. # 使用与源代码一致的BallTree参数
  70. self.normalDataBallTree = BallTree(
  71. self.matrixD,
  72. leaf_size=4,
  73. metric=lambda i, j: 1 - self.calcSimilarity(i, j) # 自定义相似度
  74. )
  75. # 使用所有数据计算残差
  76. self.healthyResidual = self.calcResidualByLocallyWeightedLR(trainDataset)
  77. return 0
  78. def calcResidualByLocallyWeightedLR(self, newStates):
  79. """优化残差计算"""
  80. if len(newStates.shape) == 1:
  81. newStates = newStates.reshape(-1, 1)
  82. dist, iList = self.normalDataBallTree.query(
  83. newStates,
  84. k=min(10, len(self.matrixD)),
  85. return_distance=True
  86. )
  87. weights = 1 / (dist + 1e-5)
  88. weights /= weights.sum(axis=1)[:, np.newaxis]
  89. est_X = np.sum(weights[:, :, np.newaxis] * self.matrixD[iList[0]], axis=1)
  90. return est_X - newStates
  91. def calcSPRT(self, newsStates, feature_weight, alpha=0.1, beta=0.1, decisionGroup=1):
  92. """优化SPRT计算"""
  93. try:
  94. stateResidual = self.calcResidualByLocallyWeightedLR(newsStates)
  95. weightedStateResidual = np.dot(stateResidual, feature_weight)
  96. weightedHealthyResidual = np.dot(self.healthyResidual, feature_weight)
  97. mu0 = np.mean(weightedHealthyResidual)
  98. sigma0 = np.std(weightedHealthyResidual)
  99. # 处理标准差为零的情况
  100. if sigma0 < 1e-5:
  101. sigma0 = 1.0 # 设为安
  102. # 向量化计算
  103. n = len(newsStates)
  104. if n < decisionGroup:
  105. return [50] # 中性值
  106. rolling_mean = np.convolve(weightedStateResidual, np.ones(decisionGroup) / decisionGroup, 'valid')
  107. si = (rolling_mean - mu0) * (rolling_mean + mu0 - 2 * mu0) / (2 * sigma0 ** 2)
  108. lowThres = np.log(beta / (1 - alpha))
  109. highThres = np.log((1 - beta) / alpha)
  110. si = np.clip(si, lowThres, highThres)
  111. si = np.where(si > 0, si / highThres, si / lowThres)
  112. flag = 100 - si * 100
  113. # 填充不足的部分
  114. if len(flag) < n:
  115. flag = np.pad(flag, (0, n - len(flag)), mode='edge')
  116. return flag.tolist()
  117. except Exception as e:
  118. logger.error(f"SPRT计算错误: {str(e)}")
  119. # 返回中性值
  120. return [50] * len(newsStates)
  121. def CRITIC_prepare(self, data, flag=1):
  122. """标准化处理"""
  123. data = data.astype(float)
  124. numeric_cols = data.select_dtypes(include=[np.number]).columns
  125. # 处理全零或常数列
  126. for col in numeric_cols:
  127. # 所有值相同
  128. if data[col].nunique() == 1:
  129. # 设为中性值
  130. data[col] = 0.5
  131. continue
  132. # 负向标准化(温度等指标)
  133. negative_cols = [col for col in numeric_cols if 'temperature' in col]
  134. for col in negative_cols:
  135. col_min = data[col].min()
  136. col_max = data[col].max()
  137. range_val = col_max - col_min
  138. if range_val < 1e-5: # 防止除零
  139. range_val = 1.0
  140. data[col] = (col_max - data[col]) / range_val
  141. # 正向标准化(其他指标)
  142. positive_cols = list(set(numeric_cols) - set(negative_cols))
  143. for col in positive_cols:
  144. col_min = data[col].min()
  145. col_max = data[col].max()
  146. range_val = col_max - col_min
  147. # 防止除零
  148. if range_val < 1e-5:
  149. range_val = 1.0
  150. data[col] = (data[col] - col_min) / range_val
  151. return data
  152. def CRITIC(self, data):
  153. """CRITIC权重计算"""
  154. data_norm = self.CRITIC_prepare(data.copy())
  155. std = data_norm.std(ddof=0).clip(lower=0.01)
  156. corr = np.abs(np.corrcoef(data_norm.T))
  157. np.fill_diagonal(corr, 0)
  158. conflict = np.sum(1 - corr, axis=1)
  159. info = std * conflict
  160. weights = info / info.sum()
  161. return pd.Series(weights, index=data.columns)
  162. def ahp(self, matrix):
  163. """AHP权重计算"""
  164. eigenvalue, eigenvector = np.linalg.eig(matrix)
  165. max_idx = np.argmax(eigenvalue)
  166. weight = eigenvector[:, max_idx].real
  167. return weight / weight.sum(), eigenvalue[max_idx].real
  168. return MSETCore()
  169. def assess_turbine(self, engine_code, data, mill_type, wind_turbine_name):
  170. """评估单个风机
  171. """
  172. results = {
  173. "engine_code": engine_code,
  174. "wind_turbine_name": wind_turbine_name,
  175. "mill_type": mill_type,
  176. "total_health_score": None,
  177. "subsystems": {},
  178. "assessed_subsystems": []
  179. }
  180. # 各子系统评估
  181. subsystems_to_assess = [
  182. ('generator', self.subsystem_config['generator'][mill_type], 1),
  183. ('nacelle', self.subsystem_config['nacelle'], 1),
  184. ('grid', self.subsystem_config['grid'], 1),
  185. ('drive_train', self.subsystem_config['drive_train'] if mill_type == 'dfig' else None, 1)
  186. ]
  187. for subsystem, config, min_features in subsystems_to_assess:
  188. if config is None:
  189. continue
  190. features = self._get_subsystem_features(config, data)
  191. # 功能1:无论特征数量是否足够都输出结果
  192. if len(features) >= min_features:
  193. assessment = self._assess_subsystem(data[features])
  194. else:
  195. assessment = {
  196. 'health_score': -1, # 特征不足时输出'-'
  197. 'weights': {},
  198. 'message': f'Insufficient features (required {min_features}, got {len(features)})'
  199. }
  200. # 功能3:删除features内容
  201. if 'features' in assessment:
  202. del assessment['features']
  203. # 最终清理:确保没有NaN值
  204. for sys, result in results["subsystems"].items():
  205. if isinstance(result['health_score'], float) and np.isnan(result['health_score']):
  206. result['health_score'] = -1
  207. result['message'] = (result.get('message') or '') + '; NaN detected'
  208. if isinstance(results["total_health_score"], float) and np.isnan(results["total_health_score"]):
  209. results["total_health_score"] = -1
  210. results["subsystems"][subsystem] = assessment
  211. # 计算整机健康度(使用新字段名)
  212. if results["subsystems"]:
  213. # 只计算健康值为数字的子系统
  214. valid_subsystems = [
  215. k for k, v in results["subsystems"].items()
  216. if isinstance(v['health_score'], (int, float)) and v['health_score'] >= 0
  217. ]
  218. if valid_subsystems:
  219. weights = self._get_subsystem_weights(valid_subsystems)
  220. health_scores = [results["subsystems"][sys]['health_score'] for sys in valid_subsystems]
  221. results["total_health_score"] = float(np.dot(health_scores, weights))
  222. results["assessed_subsystems"] = valid_subsystems
  223. logger.info(f"评估结果:{results}")
  224. return results
  225. def _get_all_possible_features(self, assessor, mill_type, available_columns):
  226. """
  227. 获取所有可能的特征列(基于实际存在的列)
  228. 参数:
  229. assessor: HealthAssessor实例
  230. mill_type: 机型类型
  231. available_columns: 数据库实际存在的列名列表
  232. """
  233. features = []
  234. available_columns_lower = [col.lower() for col in available_columns] # 不区分大小写匹配
  235. for subsys_name, subsys_config in assessor.subsystem_config.items():
  236. # 处理子系统配置
  237. if subsys_name == 'generator':
  238. config = subsys_config.get(mill_type, {})
  239. elif subsys_name == 'drive_train' and mill_type != 'dfig':
  240. continue
  241. else:
  242. config = subsys_config
  243. # 处理固定特征
  244. if 'fixed' in config:
  245. for f in config['fixed']:
  246. if f in available_columns:
  247. features.append(f)
  248. # 处理关键词特征
  249. if 'keywords' in config:
  250. for rule in config['keywords']:
  251. matched = []
  252. include_kws = [kw.lower() for kw in rule['include']]
  253. exclude_kws = [ex.lower() for ex in rule.get('exclude', [])]
  254. for col in available_columns:
  255. col_lower = col.lower()
  256. # 检查包含关键词
  257. include_ok = all(kw in col_lower for kw in include_kws)
  258. # 检查排除关键词
  259. exclude_ok = not any(ex in col_lower for ex in exclude_kws)
  260. if include_ok and exclude_ok:
  261. matched.append(col)
  262. if len(matched) >= rule.get('min_count', 1):
  263. features.extend(matched)
  264. return list(set(features)) # 去重
  265. def _get_subsystem_features(self, config: Dict, data: pd.DataFrame) -> List[str]:
  266. """最终版特征获取方法"""
  267. available_features = []
  268. # 固定特征检查(要求至少10%非空)
  269. if 'fixed' in config:
  270. for f in config['fixed']:
  271. if f in data.columns and data[f].notna().mean() > 0.1:
  272. available_features.append(f)
  273. # logger.info(f"匹配到的固定特征: {available_features}")
  274. # 关键词特征检查
  275. if 'keywords' in config:
  276. for rule in config['keywords']:
  277. matched = [
  278. col for col in data.columns
  279. if all(kw.lower() in col.lower() for kw in rule['include'])
  280. and not any(ex.lower() in col.lower() for ex in rule.get('exclude', []))
  281. and data[col].notna().mean() > 0.1 # 数据有效性检查
  282. ]
  283. if len(matched) >= rule.get('min_count', 1):
  284. available_features.extend(matched)
  285. # logger.info(f"匹配到的关键词特征: {available_features}")
  286. return list(set(available_features))
  287. def _assess_subsystem(self, data: pd.DataFrame) -> Dict:
  288. """评估子系统(与源代码逻辑完全一致)"""
  289. # 数据清洗
  290. clean_data = data.dropna()
  291. if len(clean_data) < 20: # 数据量不足
  292. return {'health_score': -1, 'weights': {}, 'features': list(data.columns), 'message': 'Insufficient data'}
  293. try:
  294. # 标准化
  295. normalized_data = self.mset.CRITIC_prepare(clean_data)
  296. # 计算权重
  297. weights = self.mset.CRITIC(normalized_data)
  298. # MSET评估
  299. health_score = self._run_mset_assessment(normalized_data.values, weights.values)
  300. return {
  301. 'health_score': float(health_score),
  302. 'weights': weights.to_dict(),
  303. 'features': list(data.columns)
  304. }
  305. except Exception as e:
  306. return {'health_score': -1, 'weights': {}, 'features': list(data.columns), 'message': str(e)}
  307. @lru_cache(maxsize=10)
  308. def _get_mset_model(self, train_data: tuple):
  309. """缓存MSET模型"""
  310. # 注意:由于lru_cache需要可哈希参数,这里使用元组
  311. arr = np.array(train_data)
  312. model = self._create_mset_core()
  313. model.genDLMatrix(arr)
  314. return model
  315. def _run_mset_assessment(self, data: np.ndarray, weights: np.ndarray) -> float:
  316. """执行MSET评估"""
  317. # 检查权重有效性
  318. if np.isnan(weights).any() or np.isinf(weights).any():
  319. # 重置为等权重
  320. weights = np.ones_like(weights) / len(weights)
  321. # 分割训练集和测试集
  322. split_idx = len(data) // 2
  323. train_data = data[:split_idx]
  324. test_data = data[split_idx:]
  325. # 使用缓存模型
  326. try:
  327. model = self._get_mset_model(tuple(map(tuple, train_data)))
  328. flags = model.calcSPRT(test_data, weights)
  329. # 过滤NaN值并计算均值
  330. valid_flags = [x for x in flags if not np.isnan(x)]
  331. if not valid_flags:
  332. # 默认中性值
  333. return 50.0
  334. return float(np.mean(valid_flags))
  335. except Exception as e:
  336. logger.error(f"MSET评估失败: {str(e)}")
  337. # 默认中性值
  338. return 50.0
  339. def _get_subsystem_weights(self, subsystems: List[str]) -> np.ndarray:
  340. """生成等权重的子系统权重向量"""
  341. n = len(subsystems)
  342. if n == 0:
  343. return np.array([])
  344. # 直接返回等权重向量
  345. return np.ones(n) / n