MSET_Temp.py 11 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281
  1. import math
  2. import os
  3. import joblib
  4. import numpy as np
  5. import pandas as pd
  6. from sqlalchemy import text
  7. from sklearn.neighbors import BallTree
  8. from app.config import dataBase
  9. from app.database import get_engine
  10. from typing import Dict
  11. class MSET_Temp:
  12. """
  13. MSET + SPRT 温度分析类:
  14. - 离线训练:genDLMatrix → save_model
  15. - 在线推理:load_model → predict_SPRT
  16. """
  17. def __init__(self,
  18. windCode: str,
  19. windTurbineNumberList: list[str],
  20. startTime: str,
  21. endTime: str):
  22. self.windCode = windCode.strip()
  23. self.windTurbineNumberList = windTurbineNumberList or []
  24. self.startTime = startTime
  25. self.endTime = endTime
  26. # 离线训练/加载后赋值
  27. self.matrixD = None
  28. self.healthyResidual = None
  29. self.normalDataBallTree = None
  30. # SPRT 参数(离线训练时设置)
  31. self.feature_weight: np.ndarray | None = None
  32. self.alpha: float = 0.1
  33. self.beta: float = 0.1
  34. def _get_data_by_filter(self) -> pd.DataFrame:
  35. """
  36. 在线推理专用:根据 self.windTurbineNumberList & 时间拉数据;
  37. 如果列表为空,则拉全场数据。
  38. """
  39. table = f"{self.windCode}_minute"
  40. engine = get_engine(dataBase.DATA_DB)
  41. if self.windTurbineNumberList:
  42. turbines = ",".join(f"'{t}'" for t in self.windTurbineNumberList)
  43. cond = f"wind_turbine_number IN ({turbines}) AND time_stamp BETWEEN :start AND :end"
  44. else:
  45. cond = "time_stamp BETWEEN :start AND :end"
  46. sql = text(f""" SELECT * FROM {table} WHERE {cond} ORDER BY time_stamp ASC """)
  47. return pd.read_sql(sql, engine, params={"start": self.startTime, "end": self.endTime})
  48. def calcSimilarity(self, x: np.ndarray, y: np.ndarray, m: str = 'euc') -> float:
  49. if len(x) != len(y):
  50. return 0.0
  51. if m == 'cbd':
  52. return float(np.mean([1.0 / (1.0 + abs(p - q)) for p, q in zip(x, y)]))
  53. diffsq = np.sum((x - y) ** 2)
  54. return float(1.0 / (1.0 + math.sqrt(diffsq)))
  55. def genDLMatrix(self, trainDataset: np.ndarray, dataSize4D=100, dataSize4L=50) -> int:
  56. """
  57. 离线训练:构造 matrixD/matrixL/healthyResidual/BallTree
  58. """
  59. m, n = trainDataset.shape
  60. if m < dataSize4D + dataSize4L:
  61. return -1
  62. # Step1:每维最小/最大入 D
  63. D_idx, D = [], []
  64. for i in range(n):
  65. col = trainDataset[:, i]
  66. for idx in (np.argmin(col), np.argmax(col)):
  67. D.append(trainDataset[idx].tolist())
  68. D_idx.append(idx)
  69. # Step2:挑样本至 dataSize4D
  70. while len(D_idx) < dataSize4D:
  71. free = list(set(range(m)) - set(D_idx))
  72. scores = [(np.mean([1 - self.calcSimilarity(trainDataset[i], d) for d in D]), i)
  73. for i in free]
  74. _, pick = max(scores)
  75. D.append(trainDataset[pick].tolist())
  76. D_idx.append(pick)
  77. self.matrixD = np.array(D)
  78. # BallTree + healthyResidual
  79. self.normalDataBallTree = BallTree(
  80. self.matrixD,
  81. leaf_size=4,
  82. metric=lambda a, b: 1.0 - self.calcSimilarity(a, b)
  83. )
  84. # healthyResidual
  85. ests = []
  86. for x in trainDataset:
  87. dist, idxs = self.normalDataBallTree.query([x], k=20, return_distance=True)
  88. w = 1.0 / (dist[0] + 1e-1)
  89. w /= w.sum()
  90. ests.append(np.sum([wi * self.matrixD[j] for wi, j in zip(w, idxs[0])], axis=0))
  91. self.healthyResidual = np.array(ests) - trainDataset
  92. return 0
  93. def calcSPRT(self,
  94. newsStates: np.ndarray,
  95. feature_weight: np.ndarray,
  96. alpha: float = 0.1,
  97. beta: float = 0.1,
  98. decisionGroup: int = 5) -> list[float]:
  99. """
  100. Wald-SPRT 得分
  101. """
  102. # 新状态残差
  103. ests = []
  104. for x in newsStates:
  105. dist, idxs = self.normalDataBallTree.query([x], k=20, return_distance=True)
  106. w = 1.0 / (dist[0] + 1e-1);
  107. w /= w.sum()
  108. ests.append(np.sum([wi * self.matrixD[j] for wi, j in zip(w, idxs[0])], axis=0))
  109. resN = np.array(ests) - newsStates
  110. # 加权
  111. wN = [np.dot(r, feature_weight) for r in resN]
  112. wH = [np.dot(r, feature_weight) for r in self.healthyResidual]
  113. mu0, sigma0 = np.mean(wH), np.std(wH)
  114. low = math.log(beta / (1 - alpha));
  115. high = math.log((1 - beta) / alpha)
  116. flags = []
  117. for i in range(len(wN) - decisionGroup + 1):
  118. seg = wN[i:i + decisionGroup];
  119. mu1 = np.mean(seg)
  120. si = (sum(seg) * (mu1 - mu0) / sigma0 ** 2
  121. - decisionGroup * ((mu1 ** 2 - mu0 ** 2) / (2 * sigma0 ** 2)))
  122. si = max(min(si, high), low)
  123. flags.append(si / high if si > 0 else si / low)
  124. return flags
  125. def predict_SPRT(self,
  126. newsStates: np.ndarray,
  127. decisionGroup: int = 5) -> list[float]:
  128. """
  129. 在线推理:用离线保存的 matrixD/healthyResidual/feature_weight/alpha/beta
  130. """
  131. return self.calcSPRT(
  132. newsStates,
  133. self.feature_weight,
  134. alpha=self.alpha,
  135. beta=self.beta,
  136. decisionGroup=decisionGroup
  137. )
  138. def save_model(self, path: str):
  139. """
  140. Save matrixD, healthyResidual, feature_weight, alpha, beta
  141. """
  142. os.makedirs(os.path.dirname(path), exist_ok=True)
  143. joblib.dump({
  144. 'matrixD': self.matrixD,
  145. 'healthyResidual': self.healthyResidual,
  146. 'feature_weight': self.feature_weight,
  147. 'alpha': self.alpha,
  148. 'beta': self.beta,
  149. }, path)
  150. @classmethod
  151. def load_model(cls, path: str) -> 'MSET_Temp':
  152. """
  153. Load + rebuild BallTree
  154. """
  155. data = joblib.load(path)
  156. inst = cls('', [], '', '')
  157. inst.matrixD = data['matrixD']
  158. inst.healthyResidual = data['healthyResidual']
  159. inst.feature_weight = data['feature_weight']
  160. inst.alpha = data['alpha']
  161. inst.beta = data['beta']
  162. inst.normalDataBallTree = BallTree(
  163. inst.matrixD,
  164. leaf_size=4,
  165. metric=lambda a, b: 1.0 - inst.calcSimilarity(a, b)
  166. )
  167. return inst
  168. def query_surrounding_data(self, timestamp: str, minutes_around: int = 250) -> Dict:
  169. """
  170. 查询指定时间点前后50个点的数据
  171. 参数:
  172. timestamp: 中心时间点,格式为 'yyyy-mm-dd HH:MM:SS'
  173. minutes_around: 查询前后多少分钟的数据
  174. 返回:
  175. {
  176. 'record_count': int,
  177. 'records': List[Dict],
  178. 'columns_mapping': Dict[str, str] # 字段中英文映射
  179. }
  180. """
  181. # 中英文映射字典
  182. cn_map = {
  183. 'wind_turbine_name':'风机名称',
  184. 'time_stamp': '时间',
  185. 'active_power': '有功功率(kW)',
  186. 'rotor_speed': '风轮转速(rpm)',
  187. 'generator_speed':'发电机转速(rpm)',
  188. 'wind_velocity': '风速(m/s)',
  189. 'pitch_angle_blade_1':'桨距角1(°)',
  190. 'pitch_angle_blade_2':'桨距角2(°)',
  191. 'pitch_angle_blade_3':'桨距角3(°)',
  192. 'cabin_position':'机舱位置(°)',
  193. 'true_wind_direction':'绝对风向(°)',
  194. 'yaw_error1':'对风角度(°)',
  195. 'set_value_of_active_power':'有功功率设定值(kW)',
  196. 'gearbox_oil_temperature':'齿轮箱油温(℃)',
  197. 'generatordrive_end_bearing_temperature':'发电机驱动端轴承温度(℃)',
  198. 'generatornon_drive_end_bearing_temperature':'发电机非驱动端轴承温度(℃)',
  199. 'cabin_temperature':'机舱内温度(℃)',
  200. 'twisted_cable_angle':'扭缆角度(°)',
  201. 'outside_cabin_temperature':'环境温度(℃)',
  202. 'main_bearing_temperature':'主轴承轴承温度(℃)',
  203. 'main_bearing_temperature_2': '主轴承轴承温度2(℃)',
  204. 'gearbox_high_speed_shaft_bearing_temperature':'齿轮箱高速轴轴承温度(℃)',
  205. 'gearboxmedium_speed_shaftbearing_temperature':'齿轮箱中速轴轴承温度(℃)',
  206. 'gearbox_low_speed_shaft_bearing_temperature':'齿轮箱低速轴轴承温度(℃)',
  207. 'generator_winding1_temperature':'发电机绕组1温度(℃)',
  208. 'generator_winding2_temperature':'发电机绕组2温度(℃)',
  209. 'generator_winding3_temperature':'发电机绕组3温度(℃)',
  210. 'grid_a_phase_current':'电网A相电流(A)',
  211. 'grid_b_phase_current': '电网B相电流(A)',
  212. 'grid_c_phase_current': '电网C相电流(A)'
  213. }
  214. table = f"{self.windCode}_minute"
  215. engine = get_engine(dataBase.DATA_DB)
  216. # 查询数据
  217. sql = text(f"""
  218. SELECT *
  219. FROM {table}
  220. WHERE wind_turbine_number IN ({','.join([f"'{t}'" for t in self.windTurbineNumberList])})
  221. AND time_stamp BETWEEN
  222. DATE_SUB(:timestamp, INTERVAL :minutes MINUTE)
  223. AND DATE_ADD(:timestamp, INTERVAL :minutes MINUTE)
  224. ORDER BY time_stamp ASC
  225. """)
  226. df = pd.read_sql(sql, engine, params={
  227. "timestamp": timestamp,
  228. "minutes": minutes_around
  229. })
  230. # 打印查询到的数据条数
  231. record_count = len(df)
  232. print(f"查询到 {record_count} 条数据")
  233. if df.empty:
  234. return {
  235. 'record_count': 0,
  236. 'records': [],
  237. 'columns_mapping': {}
  238. }
  239. # 删除空列和不需要的列
  240. cols_to_drop = ['wind_turbine_number', 'reactive_power','lab', 'year', 'month','day','year_month','front_back_vibration_of_the_cabin','side_to_side_vibration_of_the_cabin',
  241. 'actual_torque','given_torque','clockwise_yaw_count','counterclockwise_yaw_count','unusable','power_curve_available','required_gearbox_speed','inverter_speed_master_control',
  242. 'wind_turbine_status','wind_turbine_status2','turbulence_intensity'
  243. ]
  244. cols_to_drop = [col for col in cols_to_drop if col in df.columns]
  245. df = df.drop(columns=cols_to_drop)
  246. df = df.dropna(axis=1, how='all')
  247. # 转换字段名和格式
  248. df['time_stamp'] = df['time_stamp'].astype(str)
  249. records = df.rename(columns=cn_map).to_dict('records')
  250. return {
  251. 'record_count': record_count,
  252. 'records': records
  253. }