# # Day-ahead load forecasting # # DB : MS SQL # # Program Language : Python # # kgpark@hdc-icontrols.com # April 10, 2020 # ### BEMS 데이터 수집 메카니즘 # #### 데이터 별로 수집 타입에 따라 다르지만, Raw 테이블에 적산 값으로 저장이 되고 15min 테이블에서 해당 시간대와 그 전 시간대의 차이 값을 입력한다. # #### DGW 혹은 시스템에 이상이 생겼을 때, 데이터가 들어오지 않거나 0으로 입력된다. # #### 1시간 테이블은 15분 테이블에서 각 15분, 30분, 45분, 60분의 데이터 합산 값이 나왔다. # #### 합산 값으로 저장되다보니 4개 포인트 중 적어도 하나만 있어도 1시간 데이터로 저장이 된다. # #### 따라서, 15분 데이터를 전처리하는 것이 주효하고 데이터가 없거나 0값을 검출하여 비정상 데이터로 추정하는 것을 추천한다. # #### 또한, 1시간 단위로 데이터 주기를 변환한다면 15분 테이블의 4개 포인트 중 하나라도 값을 모른다면 그 시간의 데이터가 비정상이라고 가정하는 것을 추천한다. import matplotlib.pyplot as plt import pymssql import datetime import numpy as np import math from korean_lunar_calendar import KoreanLunarCalendar import calendar import configparser import sys import time # ## Define functions ### Define day-type def getDayName(year, month, day): return ['MON','TUE','WED','THU','FRI','SAT','SUN'][datetime.date(year, month, day).weekday()] def getDayType(DateinDay, Period, SpecialHoliday): DoW=[]; # Day of Week for i in range(Period): if DateinDay[i].year==2019 and DateinDay[i].month==5 and DateinDay[i].day==18: DoW.append([5, DateinDay[i]]) elif getDayName(DateinDay[i].year,DateinDay[i].month,DateinDay[i].day) == 'MON': DoW.append([1, DateinDay[i]]) elif getDayName(DateinDay[i].year,DateinDay[i].month,DateinDay[i].day) == 'TUE': DoW.append([2, DateinDay[i]]) elif getDayName(DateinDay[i].year,DateinDay[i].month,DateinDay[i].day) == 'WED': DoW.append([3, DateinDay[i]]) elif getDayName(DateinDay[i].year,DateinDay[i].month,DateinDay[i].day) == 'THU': DoW.append([4, DateinDay[i]]) elif getDayName(DateinDay[i].year,DateinDay[i].month,DateinDay[i].day) == 'FRI': DoW.append([5, DateinDay[i]]) elif getDayName(DateinDay[i].year,DateinDay[i].month,DateinDay[i].day) == 'SAT': DoW.append([6, DateinDay[i]]) elif getDayName(DateinDay[i].year,DateinDay[i].month,DateinDay[i].day) == 'SUN': DoW.append([7, DateinDay[i]]) for j in range(len(SpecialHoliday)): if SpecialHoliday[j] == datetime.date(DateinDay[i].year,DateinDay[i].month,DateinDay[i].day): DoW[-1][0] = 8 break ### W:1, N:2, ### W: Workday, N: Non-workday DayType=[] for i in range(Period): if DoW[i][0] <= 5: DayType.append([1, DateinDay[i]]) elif DoW[i][0] > 5: DayType.append([2, DateinDay[i]]) return DoW, DayType def Reconstruction(DayType, DatainHour, mark, DataRes, isRecent): ReconstructedData=[] DayType1h=[] Day_len = len(DayType) # Rearrange data in hour unit for i in range(Day_len): if i == Day_len-1 and isRecent: Time_len = len(DatainHour) - i*DataRes else: Time_len=DataRes for j in range(Time_len): DayType1h.append([DatainHour[i*DataRes + j], DayType[i][0], datetime.datetime(DayType[i][1].year, DayType[i][1].month, DayType[i][1].day, j, 0)]) ## data, daytype, time # 비정상 데이터보다 앞선 시간의 데이터 중 DayType이 같고 시간이 같은 5개 날 데이터의 평균으로 복원함 for i in reversed(range(len(DayType1h))): AccData=[] cnt=0 if math.isnan(DayType1h[i][0]): for j in range(len(DayType1h)): if cnt > 5: break if i < j and DayType1h[j][1] == DayType1h[i][1] and DayType1h[j][2].hour == DayType1h[i][2].hour and (not math.isnan(DayType1h[j][0])): AccData.append(DayType1h[j][0]) cnt += 1 DayType1h[i][0] = np.mean(AccData) ReconstructedData.append(DayType1h[i][0]) ReconstructedData.reverse() ### Double-checking for the data which is not reconstructed, especially in front for i in range(len(DayType1h)): AccData=[] cnt=0 if math.isnan(DayType1h[i][0]): #print('Here is NaN!!',ReconstructedData[i],i,DayType1h[i][2].hour, DayType1h[i][1]) for j in reversed(range(len(DayType1h))): if cnt > 5: break if i > j and DayType1h[j][1] == DayType1h[i][1] and DayType1h[j][2].hour == DayType1h[i][2].hour and (not math.isnan(DayType1h[j][0])): AccData.append(DayType1h[j][0]) cnt += 1 ReconstructedData[i] = np.mean(AccData) return ReconstructedData, DayType1h ## For day-ahead linear prediction def lpc_pred_DayAhead(Data_trn, DayType_trn, cov_lth, DayType_tst, DataRes): # Calculating the filter bank for each hour and day-type using traing set for c_w in range(1,3): DayType_trn[0,0]=0 CP_pred_fb=np.zeros(Data_trn.shape) lpc_fb=np.zeros([cov_lth[c_w-1],DataRes]) Prv_A=[] Prv_A=np.transpose(Data_trn[DataRes-cov_lth[c_w-1]:DataRes,np.where(DayType_trn == c_w)[0]-1]) for hr_i in range(24): lpc_fb[:,hr_i]=np.dot(np.linalg.pinv(Prv_A), np.transpose(Data_trn[hr_i,np.where(DayType_trn == c_w)[0]])) if c_w == 1: lpc_fb1=lpc_fb elif c_w == 2: lpc_fb2=lpc_fb ## For testing if DayType_tst[0,0] == 1: lpc_t=lpc_fb1 elif DayType_tst[0,0] == 2: lpc_t=lpc_fb2 Data_tt=Data_trn[:,-1] # Load prediction for test day based on the filter bank CP_pred=np.transpose(np.dot(np.transpose(Data_tt[DataRes-cov_lth[DayType_tst[0,0]-1]:DataRes+1]),lpc_t)) return CP_pred ## For step-ahead linear prediction def lpc_pred_OneStepAhead(Data_trn, DayType_trn, cov_lth, DayType_tst, DataRes): for c_w in range(1,3): DayType_trn[0,0]=0 lpc_fb=np.zeros([cov_lth[c_w-1],DataRes]) Prv_A=[] Prv_A=np.transpose(Data_trn[DataRes-cov_lth[c_w-1]:DataRes,np.where(DayType_trn == c_w)[0]-1]) lpc_fb=np.dot(np.linalg.pinv(Prv_A), np.transpose(Data_trn[0,np.where(DayType_trn == c_w)[0]])) if c_w == 1: lpc_fb1=lpc_fb elif c_w == 2: lpc_fb2=lpc_fb ## Testing if DayType_tst[0,0] == 1: lpc_t=lpc_fb1 elif DayType_tst[0,0] == 2: lpc_t=lpc_fb2 Data_tt=Data_trn[:,-1] CP_pred=np.transpose(np.dot(np.transpose(Data_tt[DataRes-cov_lth[DayType_tst[0,0]-1]:DataRes+1]),lpc_t)) return CP_pred ## Measure def MAPE(y_observed, y_pred): return np.mean(np.abs((y_observed - y_pred) / y_observed)) * 100 def MAE(y_observed, y_pred): return np.mean(np.abs(y_observed - y_pred)) def MBE(y_observed, y_pred): return (np.sum((y_observed - y_pred))/(len(y_observed)*np.mean(y_observed)))*100 def CVRMSE(y_observed, y_pred): return (np.sqrt(np.mean((y_observed - y_pred)*(y_observed - y_pred)))/np.mean(y_observed))*100 ## Check for normal time stamp def Check_AlivedTimeStamp(RawData, ComparedData, idx_raw, idx_comp): if datetime.date(RawData[idx_raw][4].year,RawData[idx_raw][4].month,RawData[idx_raw][4].day) == datetime.date(ComparedData[idx_comp].year, ComparedData[idx_comp].month, ComparedData[idx_comp].day) and datetime.time(RawData[idx_raw][4].hour,RawData[idx_raw][4].minute) == datetime.time(ComparedData[idx_comp].hour, ComparedData[idx_comp].minute): isAlived = True else: isAlived = False return isAlived if __name__ == "__main__" : ## Check every hour on the hour operating infinite loop while True: now = datetime.datetime.now().now() ## distinguish real time update and specific day ## 자정에 생기는 인덱싱 문제로 0시에는 16분에 업데이트 if (now.hour != 0 and now.minute == 1) or (now.hour == 0 and now.minute == 16): PredctionActive = True else: PredctionActive = False if now.second > 55: print("[ Current Time -", now.hour,":", now.minute,":", now.second,"], " "Sleeping for 30 seconds... Prediction starts every hour") time.sleep(30) else: print("[ Current Time -", now.hour,":", now.minute,":", now.second,"], " "Sleeping for 60 seconds... Prediction starts every hour") time.sleep(60) if PredctionActive: ## Loading .ini file myINI = configparser.ConfigParser() myINI.read("Config.ini", "utf-8" ) # MSSQL Access conn = pymssql.connect(host=myINI.get('LocalDB_Info','ip_address'), user=myINI.get('LocalDB_Info','user_id'), password=myINI.get('LocalDB_Info','user_password'), database=myINI.get('LocalDB_Info','db_name'), autocommit=True) # Create Cursor from Connection cursor = conn.cursor() # Execute SQL (Electric consumption) cursor.execute('SELECT * FROM BemsConfigData where SiteId = 1') rowDB_info = cursor.fetchone() conn.close() loadDBIP = rowDB_info[1] loadDBUserID = rowDB_info[2] loadDBUserPW = rowDB_info[3] loadDBName = rowDB_info[4] targetDBIP = rowDB_info[5] targetDBUserID = rowDB_info[6] targetDBUserPW = rowDB_info[7] targetDBName = rowDB_info[8] linearFilterLength = '24,24' print("=================== Prediction start! ===================") startday = datetime.date(2019,1,1) # ## Data accumulation isRecent = True lastday = datetime.date(now.year, now.month, now.day) if startday < datetime.date(2017,1,1): print('[ERROR] 데이터 최소 시작 시점은 2017.01.01 입니다') elif startday > lastday: print('[ERROR] 예측 타깃 시작시점이 데이터 시작 시점보다 작을 수 없습니다') now_ = datetime.date(now.year, now.month, now.day) # 학습데이터의 기간은 최대 2년으로 한정 if (startday-now_).days > 730: Ago_2year = now_ + timedelta(days=-730) startday = datetime.date(Ago_2year.year, Ago_2year.month, Ago_2year.day) # MSSQL Access conn = pymssql.connect(host = loadDBIP, user = loadDBUserID, password = loadDBUserPW, database = loadDBName, autocommit=True) # Create Cursor from Connection cursor = conn.cursor() # Execute SQL (Electric consumption) cursor.execute('SELECT * FROM BemsMonitoringPointHistory15min where SiteId = 1 and FacilityTypeId = 99 and FacilityCode = 4863 and PropertyId = 1 order by CreatedDateTime desc') # 데이타 하나씩 Fetch하여 출력 row = cursor.fetchone() DataRes_org=96 DataRes_24=24 rawData=[] while row: row = cursor.fetchone() if datetime.date(row[4].year,row[4].month,row[4].day) < startday: break rawData.append(row) rawData.reverse() # 오름차순 정렬 # 연결 끊기 conn.close() # 현장 데이터가 없을 경우 예외처리 if now.hour == 0: hour_calib = 0 else: hour_calib = 1 if datetime.datetime(now.year, now.month, now.day, now.hour, 0, 0) - datetime.timedelta(hours=hour_calib) == datetime.datetime(rawData[-1][4].year, rawData[-1][4].month, rawData[-1][4].day, rawData[-1][4].hour, 0, 0): # MSSQL Access conn = pymssql.connect(host = loadDBIP, user = loadDBUserID, password = loadDBUserPW, database = loadDBName, autocommit=True) # Create Cursor from Connection cursor = conn.cursor() # SQL문 실행 (정기휴일) cursor.execute('SELECT * FROM CmHoliday where SiteId = 1 and IsUse = 1') # 데이타 하나씩 Fetch하여 출력 row = cursor.fetchone() regularHolidayData = [row] while row: row = cursor.fetchone() regularHolidayData.append(row) regularHolidayData = regularHolidayData[0:-1] # SQL문 실행 (비정기휴일) cursor.execute('SELECT * FROM CmHolidayCustom where SiteId = 1 and IsUse = 1') # 데이타 하나씩 Fetch하여 출력 row = cursor.fetchone() suddenHolidayData = [row] while row: row = cursor.fetchone() suddenHolidayData.append(row) suddenHolidayData = suddenHolidayData[0:-1] # 연결 끊기 conn.close() # 공휴일의 음력 계산 calendar_convert = KoreanLunarCalendar() SpecialHoliday = [] for i in range(lastday.year-startday.year+1): for j in range(len(regularHolidayData)): if regularHolidayData[j][3] == 1: if regularHolidayData[j][1] == 12 and regularHolidayData[j][2] == 30: ## 설 하루 전 연휴 계산을 위함 calendar_convert.setLunarDate(startday.year+i-1, regularHolidayData[j][1], regularHolidayData[j][2], False) SpecialHoliday.append(datetime.date(int(calendar_convert.SolarIsoFormat().split(' ')[0].split('-')[0]), int(calendar_convert.SolarIsoFormat().split(' ')[0].split('-')[1]), int(calendar_convert.SolarIsoFormat().split(' ')[0].split('-')[2]))) else: calendar_convert.setLunarDate(startday.year+i, regularHolidayData[j][1], regularHolidayData[j][2], False) SpecialHoliday.append(datetime.date(int(calendar_convert.SolarIsoFormat().split(' ')[0].split('-')[0]), int(calendar_convert.SolarIsoFormat().split(' ')[0].split('-')[1]), int(calendar_convert.SolarIsoFormat().split(' ')[0].split('-')[2]))) else: SpecialHoliday.append(datetime.date(startday.year+i,regularHolidayData[j][1],regularHolidayData[j][2])) for i in range(len(suddenHolidayData)): if suddenHolidayData[i][1].year >= startday.year: SpecialHoliday.append(datetime.date(suddenHolidayData[i][1].year, suddenHolidayData[i][1].month, suddenHolidayData[i][1].day)) SpecialHoliday=list(set(SpecialHoliday)) DayPeriod = (lastday - startday).days + 1 print('First day:',startday,',', 'Last Day:', lastday,',','Current Time:', now) print('Day period :', DayPeriod) # ## Find unkown/zero data (Bad data) StartTime = datetime.datetime(int(startday.strftime('%Y')), int(startday.strftime('%m')), int(startday.strftime('%d')), 0, 0, 0) TimeStamp_DayUnit = [] StandardTimeStamp = [] # Create normal time stamp for idx_day in range(DayPeriod): TimeStamp_DayUnit.append(startday + datetime.timedelta(days=idx_day)) if isRecent and idx_day == DayPeriod-1: if now.hour == 0: # 예외처리용 (자정에 Day count가 안되는 현상) tmp_len = 1 else: tmp_len = now.hour*4 + int(now.minute/15) for idx_time in range(tmp_len): StandardTimeStamp.append(StartTime) StartTime += datetime.timedelta(minutes = 15) else: for idx_time in range(DataRes_org): StandardTimeStamp.append(StartTime) StartTime += datetime.timedelta(minutes = 15) RawDate=[] # raw data (date) RawElectricLoad=[] # raw data (electric load) for i in range(len(rawData)): if datetime.date(rawData[i][4].year,rawData[i][4].month,rawData[i][4].day) >= startday: if datetime.date(rawData[i][4].year,rawData[i][4].month,rawData[i][4].day) <= lastday: RawDate.append(rawData[i][4]) RawElectricLoad.append(rawData[i][5]) if datetime.date(rawData[i][4].year,rawData[i][4].month,rawData[i][4].day) > lastday: break Data_len=len(RawDate) if isRecent: DataAct_len = (DayPeriod-1)*DataRes_org + now.hour*4 + int(now.minute/15) else: DataAct_len = DayPeriod*DataRes_org ### Unknown/zero data counts DataCount=[] for i in range(len(TimeStamp_DayUnit)): cnt_unk=0 # For Unknown data count cnt_zero=0 # zero data count for j in range(Data_len): if TimeStamp_DayUnit[i] == datetime.date(RawDate[j].year,RawDate[j].month,RawDate[j].day): cnt_unk += 1 if RawElectricLoad[j] == 0: cnt_zero += 1 if isRecent and i==len(TimeStamp_DayUnit)-1: DataCount.append([TimeStamp_DayUnit[i], now.hour*4 + int(now.minute/15) - cnt_unk, cnt_zero]) else: DataCount.append([TimeStamp_DayUnit[i], DataRes_org-cnt_unk, cnt_zero]) ## Visualization ## 월 인덱스 설정 ## idxCal=[] idxCalName=[] idxCal.append(0) for y_idx in range(lastday.year - startday.year + 1): if startday.year == lastday.year: for m_idx in range(lastday.month - startday.month + 1): month = startday.month + m_idx idxCal.append(idxCal[-1] + calendar.monthrange(startday.year, month)[1]) idxCalName.append(calendar.month_name[month]) else: if y_idx == 0: ## 첫번째 해 for m_idx in range(13-startday.month): month = startday.month + m_idx idxCal.append(idxCal[-1] + calendar.monthrange(startday.year, month)[1]) idxCalName.append(calendar.month_name[month]) elif y_idx !=0 and y_idx == lastday.year - startday.year: ## 마지막 해 for m_idx in range(lastday.month): month = m_idx + 1 idxCal.append(idxCal[-1] + calendar.monthrange(lastday.year, month)[1]) idxCalName.append(calendar.month_name[month]) else: for m_idx in range(12): month = m_idx + 1 idxCal.append(idxCal[-1] + calendar.monthrange(startday.year+y_idx, month)[1]) idxCalName.append(calendar.month_name[month]) DataCountMat=np.matrix(DataCount) print("The number of unknown data:",sum(DataCountMat[:,1]), ", The number of zero data:", sum(DataCountMat[:,2])) plt.figure(figsize=(16,9)) plt.subplot(311) plt.plot(DataCountMat[:,1],label='Unknown data', linewidth = 2) plt.plot(DataCountMat[:,2],label='Zero data', linewidth = 2) # plt.xlabel('Months', fontsize = 16) plt.ylabel('Data counts', fontsize = 14) plt.legend(loc='upper left', fontsize = 14) plt.title("Unknown/zero electric load data per 15min. unit ("+str(startday.year)+"."+str(startday.month)+"."+str(startday.day)+" - "+str(lastday.year)+"."+str(lastday.month)+"."+str(lastday.day)+")", fontsize = 14) plt.xlim(idxCal[0], idxCal[-1]) plt.xticks(idxCal, idxCalName, fontsize=6.5) plt.yticks(fontsize=14) print("Bad data detection complete!") ### NaN-padding after finding unknown data ######## 현재 DB 특성상 값이 0으로 찍히거나 시간테이블의 행 자체가 없는 경우가 있고, 이 데이터 1 step 앞뒤로 데이터가 비정상일 확률이 높으므로 비정상데이터 뿐만 아니라 앞뒤 1 step까지 nan으로 처리함 ElectricLoad_Un_ZP=[] RawDate=[] idx=0 idx2=0 isBadData = False for i in range(DataAct_len): if datetime.date(rawData[idx][4].year,rawData[idx][4].month,rawData[idx][4].day) >= startday and datetime.date(rawData[idx][4].year,rawData[idx][4].month,rawData[idx][4].day) <= lastday: RawDate.append(StandardTimeStamp[idx2]) if isBadData == True: ElectricLoad_Un_ZP.append(np.nan) isBadData=False elif rawData[idx][5]==0: ElectricLoad_Un_ZP[-1]=np.nan ElectricLoad_Un_ZP.append(np.nan) if rawData[idx+1][5] > 0 and Check_AlivedTimeStamp(rawData, StandardTimeStamp, idx+1, idx2+1): isBadData = True elif Check_AlivedTimeStamp(rawData, StandardTimeStamp, idx, idx2): ElectricLoad_Un_ZP.append(rawData[idx][5]) else: ElectricLoad_Un_ZP[-1]=np.nan ElectricLoad_Un_ZP.append(np.nan) if rawData[idx+1][5] > 0 and Check_AlivedTimeStamp(rawData, StandardTimeStamp, idx+1, idx2+1): isBadData = True idx -= 1 idx2 += 1 idx += 1 print('NaN-padding complete!') # ## Decimation to 1-hour period ElectricLoad_1h = [] for i in range(DayPeriod): if i == DayPeriod-1 and isRecent: Time_len = DataAct_len - i*DataRes_org + 1 else: Time_len = DataRes_org isNaN=False for j in range(Time_len): if ElectricLoad_Un_ZP[i*4 + j] == np.nan: isNaN=True if j%4==3: if isNaN: ElectricLoad_1h.append(np.nan) else: ElectricLoad_1h.append(sum(ElectricLoad_Un_ZP[i*DataRes_org + j-3:i*DataRes_org + j+1])) print('Decimation to 1hour complete!') # ## Data reconstruction using similar-day approach DateinDay=[] for k in range(DayPeriod): DateinDay.append(RawDate[k*DataRes_org]) DoW, DayType = getDayType(DateinDay, DayPeriod, SpecialHoliday) # Find the similar-day and reconstructed data marking=np.nan ReconstructedData, DayType1h = Reconstruction(DayType, ElectricLoad_1h, marking, DataRes_24, isRecent) plt.subplot(312) plt.plot(ReconstructedData, '*-', label='Reconstructed data',linewidth=3) plt.plot(ElectricLoad_1h, '--', label='Raw data',linewidth=3) plt.legend(loc='upper right', fontsize = 14) plt.ylabel('Power [kW]', fontsize = 14) plt.yticks(fontsize=14) plt.xticks([0],fontsize=14) plt.xlim((DayPeriod-10)*24, DayPeriod*24) plt.title('Raw & reconstructed data in the latest 10 days',fontsize=14) print('Reconstruct complete!') # ## Day-ahead load forecasting ####### Convert to matrix ReconstructedData_Arr=np.zeros((DataRes_24, DayPeriod)) for i in range(DayPeriod): if isRecent and i==DayPeriod-1: for j in range(len(ReconstructedData)%DataRes_24): ReconstructedData_Arr[j,i]=ReconstructedData[i*DataRes_24+j] else: for j in range(DataRes_24): ReconstructedData_Arr[j,i]=ReconstructedData[i*DataRes_24+j] trn_period=DayPeriod - 1 DayType_m=np.matrix(DayType) Data_trn=ReconstructedData_Arr[:,0:trn_period] Data_tst=ReconstructedData_Arr[:,trn_period] DayType_trn=DayType_m[0:trn_period,:] DayType_tst=DayType_m[trn_period,:] cov_lth=np.array([int(linearFilterLength.split(',')[0]),int(linearFilterLength.split(',')[1])]) y_pred_dayAhead = lpc_pred_DayAhead(Data_trn, DayType_trn, cov_lth, DayType_tst, DataRes_24) print('-------------------------Day-ahead prediction result-------------------------') if isRecent: if now.hour == 0: print('MAPE :', MAPE(Data_tst[0],y_pred_dayAhead[0]), 'MAE :', MAE(Data_tst[0],y_pred_dayAhead[0])) else: print('MAPE :', MAPE(Data_tst[0:now.hour],y_pred_dayAhead[0:now.hour]), 'MAE :', MAE(Data_tst[0:now.hour],y_pred_dayAhead[0:now.hour])) else: print('MAPE :', MAPE(Data_tst,y_pred_dayAhead),'MAE :', MAE(Data_tst,y_pred_dayAhead)) print('MBE :', MBE(Data_tst,y_pred_dayAhead), 'CVRMSE :', CVRMSE(Data_tst,y_pred_dayAhead)) print('-------------------------------------------------------------------------------') # ## One-step-ahead load forecasting y_pred_oneStep=[] Data_tst_oneStep=[] if isRecent: dayHour = now.hour + 1 else: dayHour = DataRes_24 for i in range(dayHour): ####### Convert to matrix ReconstructedData_tmp=ReconstructedData[i:] if isRecent: for ii in range(DataRes_24-i): ReconstructedData_tmp.append(np.nan) for ii in range(i): ReconstructedData_tmp.append(np.nan) ReconstructedData_Arr_oneStep=np.zeros((DataRes_24, DayPeriod)) for j in range(DayPeriod): for k in range(DataRes_24): ReconstructedData_Arr_oneStep[k,j]=ReconstructedData_tmp[j*DataRes_24+k] Data_trn=ReconstructedData_Arr_oneStep[:,0:trn_period] if isRecent: Data_tst_oneStep.append(ReconstructedData_Arr_oneStep[i,trn_period]) else: Data_tst_oneStep=ReconstructedData_Arr[:,trn_period] y_pred_oneStep.append(lpc_pred_OneStepAhead(Data_trn, DayType_trn, cov_lth, DayType_tst, DataRes_24)) print('-------------------------OneStep-ahead prediction result-------------------------') if isRecent: if now.hour == 0: print('MAPE :', MAPE(Data_tst[0],y_pred_oneStep[0]), 'MAE :', MAE(Data_tst[0],y_pred_oneStep[0])) else: print('MAPE :', MAPE(Data_tst[0:now.hour],y_pred_oneStep[0:now.hour]), 'MAE :', MAE(Data_tst[0:now.hour],y_pred_oneStep[0:now.hour])) else: print('MAPE :', MAPE(Data_tst_oneStep,y_pred_oneStep),'MAE :', MAE(Data_tst_oneStep,y_pred_oneStep)) print('-------------------------------------------------------------------------------') plt.subplot(313) plt.grid(b=True, which='both',axis='y') if isRecent: plt.plot(ReconstructedData_Arr[0:now.hour,trn_period], label='Observed data', linewidth=3) else: plt.plot(ReconstructedData_Arr[:,trn_period], label='Observed data', linewidth=3) plt.plot(y_pred_dayAhead, '--', label='Day-ahead Prediction', linewidth=3) plt.plot(y_pred_oneStep, '*-.', label='OneStep-ahead Prediction', MarkerSize=10, linewidth=3) plt.xlabel('Time [hour]', fontsize = 14) plt.ylabel('Power [kW]', fontsize = 14) plt.legend(loc='upper right', fontsize = 14) plt.xticks([6,12,18,24],['6','12','18','24'], fontsize = 14) plt.yticks(fontsize = 14) plt.ylim(min(ReconstructedData)*0.9,max(ReconstructedData)*1.1) if isRecent: plt.title("Electric load forecasting on "+str(now.year)+"."+str(now.month)+"."+str(now.day)+" (Updated every hour) - DGB 2nd branch", fontsize=14) else: plt.title("Electric load forecasting on "+str(lastday.year)+"."+str(lastday.month)+"."+str(lastday.day)+" (Updated every hour) - DGB 2nd branch", fontsize=14) #plt.show() print("=================== Prediction was successfully finished! ===================") fig = plt.gcf() if isRecent: # Save the figure file of result # fig.savefig("Result of electric load forecasting on "+str(now.year)+"."+str(now.month)+"."+str(now.day)+" "+str(now.hour)+"h"+str(now.minute)+"m - DGB 2nd branch.png", dpi=fig.dpi) ### One-hour-ahead load forecasting updated every 1 minute # MSSQL Access conn = pymssql.connect(host = targetDBIP, user = targetDBUserID, password = targetDBUserPW, database = targetDBName, autocommit=True) # Create Cursor from Connection cursor = conn.cursor() ######################################################################### #try: # cursor.execute("INSERT INTO " + targetDBName + ".dbo.BemsMonitoringPointForecastingHourAhead (SiteId,FacilityTypeId,FacilityCode,PropertyId,CreatedDateTime,TargetDateTime,ForecastedValue) VALUES(1,99,4863,1,'" + datetime.datetime.now().strftime('%Y-%m-%d %H:%M:%S') + "','"+ datetime.datetime(now.year,now.month,now.day,now.hour,0,0).strftime('%Y-%m-%d %H:00:00') + "',"+str(y_pred_oneStep[-1])+")") #except: # print('Hour-ahead forecasted data already exists! (' + (datetime.datetime(now.year,now.month,now.day,now.hour,0,0) - datetime.timedelta(hours=1)).strftime('%Y-%m-%d %H:00:00') + ')') ######################################################################### ### Day-ahead load forecasting updated every midnight #if now.hour == 0: # # Create Cursor from Connection # cursor = conn.cursor() # for i in range(len(y_pred_dayAhead)): # try: # cursor.execute("INSERT INTO " + targetDBName + ".dbo.BemsMonitoringPointForecastingDayAhead (SiteId,FacilityTypeId,FacilityCode,PropertyId,CreatedDateTime,TargetDateTime,ForecastedValue) VALUES(1,99,4863,1,'" + datetime.datetime.now().strftime('%Y-%m-%d %H:%M:%S') + "','" + (datetime.datetime(now.year,now.month,now.day,0,0,0) + datetime.timedelta(hours=i)).strftime('%Y-%m-%d %H:00:00') + "'," + str(y_pred_dayAhead[i]) + ")") # except: # print('Day-ahead forecasted data already exists! ('+(datetime.datetime(now.year,now.month,now.day,0,0,0) + datetime.timedelta(hours=i)).strftime('%Y-%m-%d %H:00:00')+')') ######################################################################### ######################################################################### # Restore the previous data and save data # Hour-ahead (지우고 다시 쓰기) currentTime = datetime.datetime(now.year, now.month, now.day, now.hour, 0, 0) tmpTime = [datetime.datetime(now.year, now.month, now.day, 0, 0, 0)] while tmpTime[-1] < currentTime: tmpTime.append(tmpTime[-1] + datetime.timedelta(hours=1)) cursor.execute("DELETE " + targetDBName + ".dbo.BemsMonitoringPointForecastingHourAhead where SiteId = 1 and FacilityTypeId = 99 and FacilityCode = 4863 and PropertyId = 1 and TargetDateTime >= '" + datetime.datetime.now().strftime('%Y-%m-%d 00:00:00') + "' and TargetDateTime < '" + (datetime.datetime.now() + datetime.timedelta(days=1)).strftime('%Y-%m-%d 00:00:00') + "'") for i in range(len(tmpTime)): cursor.execute("INSERT INTO " + targetDBName + ".dbo.BemsMonitoringPointForecastingHourAhead (SiteId,FacilityTypeId,FacilityCode,PropertyId,CreatedDateTime,TargetDateTime,ForecastedValue) VALUES(1,99,4863,1,'" + datetime.datetime.now().strftime('%Y-%m-%d %H:%M:%S') + "','"+ datetime.datetime(tmpTime[i].year,tmpTime[i].month,tmpTime[i].day,tmpTime[i].hour,0,0).strftime('%Y-%m-%d %H:00:00') + "',"+str(y_pred_oneStep[i])+")") # Day-ahead (지우고 다시 쓰기) FinalTime = datetime.datetime(now.year, now.month, now.day, 23, 0, 0) tmpTime = [datetime.datetime(now.year, now.month, now.day, 0, 0, 0)] while tmpTime[-1] < FinalTime: tmpTime.append(tmpTime[-1] + datetime.timedelta(hours=1)) cursor.execute("DELETE " + targetDBName + ".dbo.BemsMonitoringPointForecastingDayAhead where SiteId = 1 and FacilityTypeId = 99 and FacilityCode = 4863 and PropertyId = 1 and TargetDateTime >= '" + datetime.datetime.now().strftime('%Y-%m-%d 00:00:00') + "' and TargetDateTime < '" + (datetime.datetime.now() + datetime.timedelta(days=1)).strftime('%Y-%m-%d 00:00:00') + "'") for i in range(len(tmpTime)): cursor.execute("INSERT INTO " + targetDBName + ".dbo.BemsMonitoringPointForecastingDayAhead (SiteId,FacilityTypeId,FacilityCode,PropertyId,CreatedDateTime,TargetDateTime,ForecastedValue) VALUES(1,99,4863,1,'" + datetime.datetime.now().strftime('%Y-%m-%d %H:%M:%S') + "','"+ datetime.datetime(tmpTime[i].year,tmpTime[i].month,tmpTime[i].day,tmpTime[i].hour,0,0).strftime('%Y-%m-%d %H:00:00') + "',"+str(y_pred_dayAhead[i])+")") ######################################################################### conn.close() print("The result was saved!") else: fig.savefig("Result of electric load forecasting on "+str(lastday.year)+"."+str(lastday.month)+"."+str(lastday.day)+ "- DGB 2nd branch.png", dpi=fig.dpi) plt.show() print("Sleeping for 60 seconds ...") else: print("No data ... Sleeping for 60 seconds ...") time.sleep(60)