#!/usr/bin/env python # -*- coding:utf-8 -*- # @FileName :processing_limit_power_by_wind.py # @Time :2024/12/3 14:32 # @Author :David # @Company: shenyang JY import os, time import pandas as pd import numpy as np from pymongo import MongoClient from flask import request, app from logs import Log import matplotlib.pyplot as plt import pickle, traceback current_path = os.path.dirname(__file__) parent_path = os.path.dirname(current_path) def get_data_from_mongo(args): mongodb_connection,mongodb_database,mongodb_read_table = "mongodb://root:sdhjfREWFWEF23e@192.168.1.43:30000/",args['mongodb_database'],args['mongodb_read_table'] client = MongoClient(mongodb_connection) # 选择数据库(如果数据库不存在,MongoDB 会自动创建) db = client[mongodb_database] collection = db[mongodb_read_table] # 集合名称 data_from_db = collection.find() # 这会返回一个游标(cursor) # 将游标转换为列表,并创建 pandas DataFrame df = pd.DataFrame(list(data_from_db)) client.close() return df def insert_data_into_mongo(res_df,args): mongodb_connection,mongodb_database,mongodb_write_table = "mongodb://root:sdhjfREWFWEF23e@192.168.1.43:30000/",args['mongodb_database'],args['mongodb_write_table'] client = MongoClient(mongodb_connection) db = client[mongodb_database] if mongodb_write_table in db.list_collection_names(): db[mongodb_write_table].drop() print(f"Collection '{mongodb_write_table} already exist, deleted successfully!") collection = db[mongodb_write_table] # 集合名称 # 将 DataFrame 转为字典格式 data_dict = res_df.to_dict("records") # 每一行作为一个字典 # 插入到 MongoDB collection.insert_many(data_dict) print("data inserted successfully!") @app.route('/processing_limit_power_by_wind', methods=['POST', 'GET']) def processing_limit_power_by_agcavc(): # 获取程序开始时间 start_time = time.time() result = {} success = 0 print("Program starts execution!") try: logger = Log().logger args = request.values.to_dict() weather_power = get_data_from_mongo(args) lp = LimitPower(logger, args, weather_power) weather_power = lp.clean_limited_power('') print('args', args) logger.info(args) insert_data_into_mongo(weather_power, args) success = 1 except Exception as e: my_exception = traceback.format_exc() my_exception.replace("\n", "\t") result['msg'] = my_exception end_time = time.time() result['success'] = success result['args'] = args result['start_time'] = time.strftime('%Y-%m-%d %H:%M:%S', time.localtime(start_time)) result['end_time'] = time.strftime('%Y-%m-%d %H:%M:%S', time.localtime(end_time)) print("Program execution ends!") return result class LimitPower(object): def __init__(self, logger, args, weather_power): self.logger = logger self.args = args self.weather_power = weather_power self.step = self.args.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.args.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.args.usable_power['up_offset'] down_offset = self.args.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): glob_rp = self.segment_statis() 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.args.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.args.usable_power["env"], y='C_REAL_VALUE', c='c') plt.savefig(parent_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.args.usable_power['env']]] if __name__ == '__main__': from logs import Log log = Log().logger args = {} # 实例化配置类 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()