729

itil_limited_power_cleaning.py

@@ -0,0 +1,37 @@
+import requests
+import pandas as pd
+url = 'http://itil.jiayuepowertech.com:9958/itil/api/power-limitation'
+
+
+def get_station_info(timeBegin, timeEnd):
+    url_f = url + '?timeBegin=' + timeBegin + '&timeEnd=' + timeEnd
+    res = requests.get(url_f).json()
+    ele_info = list(filter(lambda x: x['stationCode'] == 'J00557', res['data']))
+    return ele_info
+
+def cleaning_powers(power, ele_info):
+    power = power.copy()
+    start = len(power)
+    print("限电清洗前，有{}条".format(start))
+    for ele in ele_info:
+        begin = ele['timeBegin']
+        end = ele['timeEnd']
+        ele_limits = power[begin: end].index.to_list()
+        if len(ele_limits) > 0:
+            print("清洗时段：{}-{}".format(ele_limits[0], ele_limits[-1]))
+            power.drop(ele_limits, inplace=True)
+    print("限电清洗后，有{}条, 一共清洗了{}条".format(len(power), start-len(power)))
+    return power
+
+if __name__ == '__main__':
+    power = pd.read_csv('./557/power.csv', parse_dates=['C_TIME'], index_col='C_TIME')
+    dq = pd.read_csv('./557/dq.csv', parse_dates=['C_TIME'])
+
+    ele_info = get_station_info('2023-01', '2023-10-31')
+    power_filter = cleaning_powers(power, ele_info)
+    power_filter.to_csv('./557/power_filter.csv')
+
+
+    power_able = pd.merge(power, dq, on='C_TIME')
+    power_able['error'] = round(power_able['ABLE'] - power_able['C_FP_VALUE'], 3)
+    power_able.to_csv('./557/power_able.csv', index=False)

limited_power_solar.py

@@ -0,0 +1,152 @@
+import pandas as pd
+import os
+import numpy as np
+np.random.seed(42)
+import matplotlib
+matplotlib.use('Agg')
+import matplotlib.pyplot as plt
+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 = pd.merge(weather, power, on='C_TIME')
+
+    def segment_statis(self):
+        """
+        对总辐射-实际功率进行分段处理，获取分度的中位点，四分位间距和斜率
+        :return: glob_rp 总辐射分段
+        """
+        segs = [x for x in range(50, 2000, 100)]    # 对辐照度以100为间隔进行分段
+        xs = [segs[i-1]+x if i>0 else 25 for i, x in enumerate([50 for _ in segs])]  # 分段的中间点
+        glob_rp = {}       # dict: key 辐照度分段中间点 value 分段内的实际功率
+        for index, row in self.weather_power.iterrows():
+            glob_ = row[self.opt.usable_power["env"]]
+            rp = row['C_REAL_VALUE']
+            for i, seg in enumerate(segs):
+                if glob_ <= seg and not (i > 0 and rp < 1):
+                    glob_rp.setdefault(xs[i], []).append(rp)
+                    break
+        for i, x in enumerate(xs):
+            rps = glob_rp.get(x)
+            if rps is None:
+                glob_rp = {k: v for k, v in glob_rp.items() if k not in xs[xs.index(x):]}
+                break
+            x_l = xs[i-1] if i > 0 else 0
+            q2_l = glob_rp[xs[i-1]][0] if i > 0 else 0
+            q1 = np.percentile(rps, self.opt.usable_power['down_fractile'])     # 实际功率下四分位点
+            q2 = np.percentile(rps, 50)  # 实际功率中位点
+            q3 = np.percentile(rps, self.opt.usable_power['up_fractile'])     # 实际功率上四分位点
+            iqr = q3 -q1    # 四分位间距
+            k1 = round(q2/x, 5)  # 整体斜率
+            k2 = round((q2-q2_l)/(x-x_l), 5)    # 趋势斜率，相对上一个中位点
+            glob_rp[x] = [q2, iqr, k1, k2]   # 更新dict
+        return glob_rp
+
+    def mapping_relation(self, glob_rp):
+        """
+        拟合分段处理后的斜率和偏移量
+        :param glob_rp: 总辐射分段
+        :return: k_final 斜率 bias 实际功率的分布宽度， glob_rp 总辐射分段
+        """
+        ks, iqrs, delete_x, tag_x = [], [], [], []   # ks所有分段斜率集合，iqrs所有分段间距集合,delete_x删除的x坐标集合
+        for x, values in glob_rp.items():
+            k1 = values[-2]
+            k2 = values[-1]
+            iqrs.append(values[-3])
+            if k1 > 0 and k2 > 0:   # 清除趋势小于等于0的斜率
+                ks.append(k1)
+                tag_x.append(x)
+            else:
+                delete_x.append(x)
+                # print("删除的斜率：", k1, k2)
+        bias = round(np.median(iqrs), 3)  # 中位点
+        # print("++++1", ks)
+        mean = np.mean(ks)  # 均值
+        std = np.std(ks)    # 标准差
+        ks = np.array(ks)
+        z_score = (ks-mean)/std # z均值
+        # print("----", z_score)
+        outliers = np.abs(z_score) > self.opt.usable_power['outliers_threshold']    # 超过阈值为离群点
+        ks = ks[~outliers]  # 消除离群点
+        delete_x1 = list(np.array(tag_x)[outliers]) # 清除大于阈值的离群点
+        k_final = round(np.mean(ks), 5)  # 对清洗后的斜率做平均
+        # print("++++2：", ks)
+        delete_x.extend(delete_x1)
+        self.logger.info("拟合可用功率，删除的斜率：" + ' '.join([str(x) for x in delete_x]))
+        glob_rp = {k: v for k, v in glob_rp.items() if k not in delete_x}   # 清洗后剩下的分段点位
+        return k_final, bias, glob_rp
+
+    def filter_unlimited_power(self, zfs, real_power, k, b):
+        """
+        预测可用功主方法
+        :param zfs: 要预测可用功率的总辐射
+        :param k: 斜率
+        :param b: 偏移量
+        :return: 预测的可用功率
+        """
+        high = k*zfs+b/2 if k*zfs+b/2 < self.opt.cap else self.opt.cap
+        low = k*zfs-b/2 if k*zfs-b/2 > 0 else 0
+        if low <= real_power <= high:
+            return True
+        else:
+            return False
+
+    def clean_limited_power(self, name, signal=False):
+        glob_rp = self.segment_statis()
+        k_final, bias, glob_rp = self.mapping_relation(glob_rp)
+        self.opt.usable_power['k'] = float(k_final)
+        self.opt.usable_power['bias'] = float(bias)
+        new_weather_power = []
+        for index, row in self.weather_power.iterrows():
+            zfs = row[self.opt.usable_power["env"]]
+            rp = row['C_REAL_VALUE']
+            s = int(row['signal']) if signal is True else 2
+            if self.filter_unlimited_power(zfs, rp, k_final, bias * self.opt.usable_power['coe']):
+                row['c'] = 'red' if s == 0 or s == 2 else 'cornflowerblue'
+                new_weather_power.append(row)
+            else:
+                row['c'] = 'blue' if s == 1 or s == 2 else 'pink'
+                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)))
+        self.args.save_args_yml(self.opt)
+        return new_weather_power.loc[:, ['C_TIME', 'C_REAL_VALUE', 'C_ABLE_VALUE']]
+
+    def clean_limited_power_by_signal(self):
+        # powers = self.weather_power.copy()
+        self.weather_power["signal"] = self.weather_power.apply(lambda x: self.signal_result(x["C_IS_RATIONING_BY_MANUAL_CONTROL"], x["C_IS_RATIONING_BY_AUTO_CONTROL"]), axis=1)
+        # powers['c'] = self.weather_power.apply(lambda x: 'blue' if bool(x["signal"]) is True else 'red', axis=1)
+        # self.weather_power.plot.scatter(x=self.opt.usable_power["env"], y='C_REAL_VALUE', c='c')
+        # plt.savefig(current_path + '/figs/限电.png')
+        # return powers
+        new_weather_power = self.weather_power[self.weather_power['signal'] == False]
+        return new_weather_power.loc[:, ['C_TIME', 'C_REAL_VALUE', 'C_ABLE_VALUE']]
+
+    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__':
+    power = pd.read_csv('2023-12-01至2023-12-23实际功率导出文件.csv', date_parser=['时间'])
+    weather = pd.read_csv('2023-12-01至2023-12-23气象站数据导出文件.csv', date_parser=['时间'])
+    weather_power = pd.merge(weather, power, on='时间')  # 联立数据
+    # glob_rp = segment_statis(weather_power)
+    # k_final, bias, glob_rp = mapping_relation(glob_rp)

limited_power_wind.py

@@ -0,0 +1,136 @@
+import os
+import pandas as pd
+import numpy as np
+np.random.seed(42)
+from cache.data_cleaning import cleaning
+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.weather_power = self.weather_power[self.weather_power[self.opt.usable_power["env"]] >= 0]
+        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])])  # 分段的中间点
+
+    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
+        for i, x in enumerate(self.xs):
+            rps = glob_rp.get(x)
+            if rps is None:
+                continue
+            mean = np.around(np.mean(rps), 3)     # 实际功率下四分位点
+            std = np.around(np.std(rps), 3)  # 实际功率中位点
+
+            glob_rp[x] = [mean, std]   # 更新dict
+        self.saveVar(os.path.dirname(current_path)+'/var/glob_rp.pickle', glob_rp)
+        return glob_rp
+
+    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, glob_rp):
+        """
+        预测可用功主方法
+        :param zfs: 要预测可用功率的总辐射
+        :param k: 斜率
+        :param b: 偏移量
+        :return: 预测的可用功率
+        """
+        coe = self.opt.usable_power['outliers_threshold']
+        seg = self.xs[np.argmax(self.segs >= ws)]
+        if seg in glob_rp:
+            mean, std = glob_rp[seg][0], glob_rp[seg][1]
+            high = mean + std*coe if mean + std*coe < self.opt.cap else self.opt.cap
+            low = mean - std*coe if mean - std*coe > 0 else 0
+            if low <= real_power <= high:
+                return True
+            else:
+                return False
+        else:
+            return True
+
+    def clean_limited_power(self, name, signal=False):
+        glob_rp = self.segment_statis()
+        new_weather_power = []
+        # fig, ax = plt.subplots()
+        for index, row in self.weather_power.iterrows():
+            zfs = row[self.opt.usable_power["env"]]
+            rp = row['C_REAL_VALUE']
+            s = int(row['signal']) if signal is True else 2
+            if self.filter_unlimited_power(zfs, rp, glob_rp):
+                row['c'] = 'red' if s == 0 or s == 2 else 'cornflowerblue'
+                new_weather_power.append(row)
+            else:
+                row['c'] = 'blue' if s == 1 or s == 2 else 'pink'
+                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)))
+        self.args.save_args_yml(self.opt)
+        return new_weather_power.loc[:, ['C_TIME', 'C_REAL_VALUE', 'C_ABLE_VALUE']]
+
+    def clean_limited_power_by_signal(self):
+        # powers = self.weather_power.copy()
+        self.weather_power["signal"] = self.weather_power.apply(lambda x: self.signal_result(x["C_IS_RATIONING_BY_MANUAL_CONTROL"], x["C_IS_RATIONING_BY_AUTO_CONTROL"]), axis=1)
+        # powers['c'] = self.weather_power.apply(lambda x: 'blue' if bool(x["signal"]) is True else 'red', axis=1)
+        # self.weather_power.plot.scatter(x=self.opt.usable_power["env"], y='C_REAL_VALUE', c='c')
+        # plt.savefig(current_path + '/figs/限电.png')
+        new_weather_power = self.weather_power[self.weather_power['signal'] == False]
+        return new_weather_power.loc[:, ['C_TIME', 'C_REAL_VALUE', 'C_ABLE_VALUE']]
+
+    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
+
+def plt_point(glob_rp, coe):
+    """
+    对最后的可用功范围进行画图
+    :param glob_rp: 总辐射分段
+    :param width: 偏移量
+    :return:
+    """
+    fig, ax = plt.subplots()
+    for x, row in glob_rp.items():
+        ax.scatter(x, row[0]-row[1]*coe, s=10, alpha=0.5)
+        ax.scatter(x, row[0]+row[1]*coe, s=10, alpha=0.5)
+        ax.scatter(x, row[0], s=12, alpha=0.75, marker='x')
+    plt.savefig('./7月份-优化1.png')
+
+
+if __name__ == '__main__':
+    power = pd.read_csv('2023-12-01至2023-12-23实际功率导出文件.csv', date_parser=['时间'])
+    weather = pd.read_csv('2023-12-01至2023-12-23气象站数据导出文件.csv', date_parser=['时间'])
+    weather_power = pd.merge(weather, power, on='时间')  # 联立数据