|
@@ -0,0 +1,194 @@
|
|
|
|
+import copy
|
|
|
|
+import os
|
|
|
|
+import pandas as pd
|
|
|
|
+import numpy as np
|
|
|
|
+import matplotlib.pyplot as plt
|
|
|
|
+import pickle
|
|
|
|
+current_path = os.path.dirname(__file__)
|
|
|
|
+
|
|
|
|
+
|
|
|
|
+class LimitPower(object):
|
|
|
|
+ def __init__(self, logger, args, weather_power):
|
|
|
|
+ self.logger = logger
|
|
|
|
+ self.args = args
|
|
|
|
+ self.opt = self.args.parse_args_and_yaml()
|
|
|
|
+ self.weather_power = weather_power
|
|
|
|
+ self.step = self.opt.usable_power['step']
|
|
|
|
+ self.segs = np.array([x * self.step for x in range(1, 50)]) # 对风速以50为间隔进行分段
|
|
|
|
+ self.xs = np.array([self.segs[i - 1] + x if i > 0 else self.step / 2 for i, x in
|
|
|
|
+ enumerate([self.step / 2 for _ in self.segs])]) # 分段的中间点
|
|
|
|
+ self.polynomial = None
|
|
|
|
+ self.width = 0
|
|
|
|
+ self.max_ws = 50
|
|
|
|
+
|
|
|
|
+ def segment_statis(self):
|
|
|
|
+ """
|
|
|
|
+ 对机头风速-实际功率进行分段处理,获取分度的中位点,四分位间距和斜率
|
|
|
|
+ :return: glob_rp 总辐射分段
|
|
|
|
+ """
|
|
|
|
+ glob_rp = {} # dict: key 辐照度分段中间点 value 分段内的实际功率
|
|
|
|
+ for index, row in self.weather_power.iterrows():
|
|
|
|
+ ws_ = row[self.opt.usable_power["env"]]
|
|
|
|
+ rp = row['C_REAL_VALUE']
|
|
|
|
+ for i, seg in enumerate(self.segs):
|
|
|
|
+ if ws_ <= seg:
|
|
|
|
+ glob_rp.setdefault(self.xs[i], []).append(rp)
|
|
|
|
+ break
|
|
|
|
+ rms = []
|
|
|
|
+ for i, x in enumerate(self.xs):
|
|
|
|
+ rps = glob_rp.get(x)
|
|
|
|
+ if rps is None:
|
|
|
|
+ continue
|
|
|
|
+ rps = np.array(rps)
|
|
|
|
+ up = self.opt.usable_power['up_fractile']
|
|
|
|
+ down = self.opt.usable_power['down_fractile']
|
|
|
|
+ offset = 0
|
|
|
|
+ while True:
|
|
|
|
+ index = i-1 if i > 0 else 0
|
|
|
|
+ while (self.xs[index] in rms or self.xs[index] not in glob_rp) and index >0:
|
|
|
|
+ index -= 1
|
|
|
|
+ x_l = self.xs[index] if index > 0 else 0
|
|
|
|
+ q2_l = glob_rp[self.xs[index]][0] if index > 0 else 0
|
|
|
|
+ down = down + offset
|
|
|
|
+ if down > up:
|
|
|
|
+ rms.append(x)
|
|
|
|
+ self.logger.info("删除的坐标点为:{}".format(x))
|
|
|
|
+ break
|
|
|
|
+ q1 = np.percentile(rps, down) # 下四分位点
|
|
|
|
+ q2 = round(np.percentile(rps, down+(up-down)/2), 3) # 中位点
|
|
|
|
+ q3 = np.percentile(rps, up) # 上四分为点
|
|
|
|
+ # q2 = np.around(np.mean(rps[(rps >= q1) & (rps <= q3)]), 3)
|
|
|
|
+ iqr = q3 - q1 # 四分位间距
|
|
|
|
+ k2 = round((q2-q2_l)/(x-x_l), 3) # 趋势斜率
|
|
|
|
+
|
|
|
|
+ # std = np.around(np.std(rps[(rps >= q1) & (rps <= q3)]), 3)
|
|
|
|
+ # mean = np.around(np.mean(rps), 3) # 实际功率均值
|
|
|
|
+ # std = np.around(np.std(rps), 3) # 实际功率标准差
|
|
|
|
+ # print("看看q2={},mean={}".format(q2, mean))
|
|
|
|
+ if k2 >= 0:
|
|
|
|
+ glob_rp[x] = [q2, iqr] # 更新dict
|
|
|
|
+ break
|
|
|
|
+ else:
|
|
|
|
+ offset += 1
|
|
|
|
+ glob_rp = {k: glob_rp[k] for k in glob_rp.keys() if k not in rms}
|
|
|
|
+ glob_rp = {k: glob_rp[k] for k in sorted(glob_rp.keys())}
|
|
|
|
+ return glob_rp
|
|
|
|
+
|
|
|
|
+ def mapping_relation(self, glob_rp):
|
|
|
|
+ degree = self.opt.usable_power['degree']
|
|
|
|
+ xs = list(glob_rp.keys())
|
|
|
|
+ ys = [y[0] for y in glob_rp.values()]
|
|
|
|
+ self.width = np.median(np.array([y[1] for y in glob_rp.values()]))
|
|
|
|
+ coefficients = np.polyfit(xs, ys, degree)
|
|
|
|
+ self.polynomial = np.poly1d(coefficients)
|
|
|
|
+ self.max_ws = max(xs)
|
|
|
|
+ # y_fit = self.polynomial(xs)
|
|
|
|
+ # plt.scatter(xs, ys, label='Data', color='red')
|
|
|
|
+ # plt.plot(xs, y_fit, label='Fitted polynomial', color='blue')
|
|
|
|
+ # plt.plot(xs, y_fit+self.width/2, label='up polynomial', color='purple')
|
|
|
|
+ # plt.plot(xs, y_fit-self.width/2, label='down polynomial', color='green')
|
|
|
|
+ # plt.legend()
|
|
|
|
+ # plt.xlabel('x')
|
|
|
|
+ # plt.ylabel('y')
|
|
|
|
+ # plt.title(f'Polynomial Fit (degree {degree})')
|
|
|
|
+ # plt.show()
|
|
|
|
+
|
|
|
|
+ def saveVar(self, path, data):
|
|
|
|
+ os.makedirs(os.path.dirname(path), exist_ok=True)
|
|
|
|
+ with open(path, 'wb') as file:
|
|
|
|
+ pickle.dump(data, file)
|
|
|
|
+
|
|
|
|
+ def filter_unlimited_power(self, ws, real_power):
|
|
|
|
+ """
|
|
|
|
+ 预测可用功主方法
|
|
|
|
+ :param zfs: 要预测可用功率的总辐射
|
|
|
|
+ :param k: 斜率
|
|
|
|
+ :param b: 偏移量
|
|
|
|
+ :return: 预测的可用功率
|
|
|
|
+ """
|
|
|
|
+ # coe = self.opt.usable_power['outliers_threshold']
|
|
|
|
+ # seg = self.xs[np.argmax(self.segs >= ws)]
|
|
|
|
+ up_offset = self.opt.usable_power['up_offset']
|
|
|
|
+ down_offset = self.opt.usable_power['down_offset']
|
|
|
|
+ high = self.polynomial(ws) + self.width/up_offset if self.polynomial(ws) + self.width/up_offset < self.opt.cap else self.opt.cap
|
|
|
|
+ low = self.polynomial(ws) - self.width/down_offset if self.polynomial(ws) - self.width/down_offset > 0 else 0
|
|
|
|
+ if low <= real_power <= high:
|
|
|
|
+ return True
|
|
|
|
+ else:
|
|
|
|
+ return False
|
|
|
|
+
|
|
|
|
+ def clean_limited_power(self, name, cluster=False):
|
|
|
|
+ if cluster is True:
|
|
|
|
+ glob_rp = self.segment_statis()
|
|
|
|
+ self.saveVar(os.path.dirname(current_path) + '/var/glob_rp.pickle', glob_rp)
|
|
|
|
+ self.mapping_relation(glob_rp)
|
|
|
|
+ new_weather_power, number = [], 0
|
|
|
|
+ # fig, ax = plt.subplots()
|
|
|
|
+ for index, row in self.weather_power.iterrows():
|
|
|
|
+ zfs = row[self.opt.usable_power["env"]]
|
|
|
|
+ rp = row['C_REAL_VALUE']
|
|
|
|
+ if zfs < 0 or rp < 0:
|
|
|
|
+ continue
|
|
|
|
+ if self.filter_unlimited_power(zfs, rp) and zfs <= self.max_ws:
|
|
|
|
+ row['c'] = 'red'
|
|
|
|
+ new_weather_power.append(row)
|
|
|
|
+ else:
|
|
|
|
+ row['c'] = 'blue'
|
|
|
|
+ new_weather_power.append(row)
|
|
|
|
+ new_weather_power = pd.concat(new_weather_power, axis=1).T
|
|
|
|
+ new_weather_power.plot.scatter(x=self.opt.usable_power["env"], y='C_REAL_VALUE', c='c')
|
|
|
|
+ plt.savefig(current_path + '/figs/测风法{}.png'.format(name))
|
|
|
|
+ new_weather_power = new_weather_power[new_weather_power['c'] == 'red']
|
|
|
|
+ number = len(new_weather_power)
|
|
|
|
+ self.logger.info("未清洗限电前,总共有:{}条数据".format(len(self.weather_power)))
|
|
|
|
+ self.logger.info("清除限电后保留的点有:" + str(number) + " 占比:" + str(round(number / len(self.weather_power), 2)))
|
|
|
|
+ return new_weather_power.loc[:, ['C_TIME', 'C_REAL_VALUE', 'C_ABLE_VALUE', self.opt.usable_power['env']]]
|
|
|
|
+
|
|
|
|
+ def clean_limited_power_by_signal(self, name):
|
|
|
|
+ weather_power1 = self.weather_power.copy()
|
|
|
|
+ weather_power1["signal"] = weather_power1.apply(
|
|
|
|
+ lambda x: self.signal_result(x["C_IS_RATIONING_BY_MANUAL_CONTROL"], x["C_IS_RATIONING_BY_AUTO_CONTROL"]),
|
|
|
|
+ axis=1)
|
|
|
|
+ weather_power1['c'] = weather_power1.apply(lambda x: 'cornflowerblue' if bool(x["signal"]) is True else 'pink',
|
|
|
|
+ axis=1)
|
|
|
|
+ weather_power1.plot.scatter(x=self.opt.usable_power["env"], y='C_REAL_VALUE', c='c')
|
|
|
|
+ plt.savefig(current_path + '/figs/信号法{}.png'.format(name))
|
|
|
|
+ weather_power1 = weather_power1[weather_power1['signal'] == False]
|
|
|
|
+ self.logger.info("信号法-未清洗限电前,总共有:{}条数据".format(len(self.weather_power)))
|
|
|
|
+ self.logger.info("信号法-清除限电后保留的点有:" + str(len(weather_power1)) + " 占比:" + str(
|
|
|
|
+ round(len(weather_power1) / len(self.weather_power), 2)))
|
|
|
|
+ return weather_power1.loc[:, ['C_TIME', 'C_REAL_VALUE', 'C_ABLE_VALUE', self.opt.usable_power['env']]]
|
|
|
|
+
|
|
|
|
+ def signal_result(self, manual, auto):
|
|
|
|
+ if int(manual) == 0:
|
|
|
|
+ if int(auto) == 0:
|
|
|
|
+ return False
|
|
|
|
+ else:
|
|
|
|
+ return True
|
|
|
|
+ else:
|
|
|
|
+ if int(auto) == 1:
|
|
|
|
+ return True
|
|
|
|
+ else:
|
|
|
|
+ return False
|
|
|
|
+
|
|
|
|
+
|
|
|
|
+if __name__ == '__main__':
|
|
|
|
+ from logs import Log
|
|
|
|
+ from config import myargparse
|
|
|
|
+ log = Log().logger
|
|
|
|
+ # 实例化配置类
|
|
|
|
+ args = myargparse(discription="场站端配置", add_help=False)
|
|
|
|
+ power = pd.read_csv('./data/power.csv')
|
|
|
|
+ weather = pd.read_csv('./data/tower-1-process.csv')
|
|
|
|
+ weather_power = pd.merge(weather, power, on='C_TIME') # 联立数据
|
|
|
|
+ lp = LimitPower(log, args, weather_power)
|
|
|
|
+ # glob_rp = lp.segment_statis()
|
|
|
|
+ # lp.mapping_relation(glob_rp)
|
|
|
|
+ lp.clean_limited_power('测试1')
|
|
|
|
+ # glob_rp = {k: glob_rp[k] for k in sorted(glob_rp.keys())}
|
|
|
|
+ # keys = list(glob_rp.keys())
|
|
|
|
+ # values = [v[0] for v in glob_rp.values()]
|
|
|
|
+ # import matplotlib.pyplot as plt
|
|
|
|
+ # fig, ax = plt.subplots()
|
|
|
|
+ # ax.plot(keys, values)
|
|
|
|
+ # plt.show()
|