ゲリラ豪雨が発生する日の太陽光発電量を予測する

太陽光発電量の予測について、以前の記事でご紹介致しました。

https://qiita.com/kccs_api1/items/360d6b135030f46eb5b1

太陽光発電量予測は、蓄電池の自動制御など電力制御のためのインプット情報として活用されていることが多いのではないかと思います。また、ピークシフトを考慮した蓄電池制御をしている場合には、30分単位での細かい太陽光発電量予測が必要となります。

そこで、ゲリラ豪雨が発生したときに太陽光発電量はどうなるのか、予測は可能なのか、を全2回の記事で見ていきたいと思います。
ここでは、ゲリラ豪雨を「長時間予報では予測できない短時間の集中的な（1時間10mm以上）の降雨量」と定義したいと思います。
今回はその第1回目です。まずは、ゲリラ豪雨が発生したときに太陽光発電量はどうなるのか、確認したいと思います。


###ゲリラ豪雨が発生した日を特定してみる
まずは、予測では雨が降らない予定だったのに、実際は雨が降った日を特定してみます。

①予測データ
降雨量の予測は、KCCS APIサービスの過去気象予報をダウンロードをします。
ダウンロード方法は、冒頭でご紹介した以前の記事を参考にしてください。

②実績データ
KCCS APIサービスでは実績は取り扱っていません。※別記事にて、別途ご説明したいと思います。
そこで、「速報版解析雨量データ配信機能（作成頻度を10分間隔とした1kmメッシュの前１時間降水量分布）」を使いたいと思います。
本機能では、過去のデータを一括して取得する機能はないため、取得したデータを蓄積しました。

####予測及び解析データの加工
#####①予測(過去気象予報)データの加工
・あらかじめ、過去気象予報データを「WEATHER_FORCAST_RESULT.csv」として配置する。
・過去気象予報データの降雨量については、1時間ごとの後1時間の積算値であるため、9:00～18:00のデータのみとするためには、8:00～17:00のデータに絞る必要がある。
・予測対象日時の重複を排除するため、今回は、当日の予測演算実施時間が0時のデータのみのデータに絞る。

```py
#過去気象予報データの読み込み
df_forecast = pd.read_csv("WEATHER_FORCAST_RESULT.csv")
df_forecast = df_forecast.sort_values(["日時","予報演算実施日時"])
df_forecast["予報演算実施日時"] = pd.to_datetime(df_forecast["予報演算実施日時"])
df_forecast["日時"] = pd.to_datetime(df_forecast["日時"])

#予測演算実施時間が0時のみのデータにしてみる
#日時が8時から18時のみのデータにしてみる
df_forecast['予報演算実施日'] = df_forecast['予報演算実施日時'].apply(lambda x : x.strftime('%Y/%m/%d'))
df_forecast['予報演算実施時'] = df_forecast['予報演算実施日時'].apply(lambda x : x.strftime('%H:%M:%S'))
df_forecast['日'] = df_forecast['日時'].apply(lambda x : x.strftime('%Y/%m/%d'))
df_forecast['時'] = df_forecast['日時'].apply(lambda x : x.strftime('%H:%M:%S'))

df_forecast = df_forecast.dropna().reset_index(drop=True)

#0時の予測演算実施時のデータに絞る
df_forecast_0hour = df_forecast[(df_forecast['予報演算実施時']=="00:00:00") & (df_forecast['予報演算実施日']==df_forecast['日'])].reset_index(drop=True)
df_forecast_0hour = df_forecast_0hour[(df_forecast_0hour['時']>="08:00:00") & (df_forecast_0hour['時']<="17:00:00")].reset_index(drop=True)
df_forecast_0hour = df_forecast_0hour[["日時", "日", "降雨量（mm/h）", "日射量（W/m^2）","中層雲量（%）"]].copy()
df_forecast_0hour = df_forecast_0hour.rename(columns={'日時':'timestamp', '日':'date',\
                                    '降雨量（mm/h）':'rain_forecast_0h','日射量（W/m^2）':'isolation_forecast_0h', '中層雲量（%）':'cloudy_forecast_0h'})
```
出力結果は以下の通りです。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1194651/f50a5256-2707-4ebf-69d7-1866031e4a8e.png)


#####②解析（速報版解析雨量）データの加工
・あらかじめ、データフレーム「df_rain_analysis_q_short」に、「速報版解析雨量データ」を格納する。
・速報版解析雨量については、10分ごとの前1時間の積算値であるため、8:00～18:00のデータのみとするためには、9:00～18:00のデータに絞る必要がある。
・10分ごとにデータが出力されるため、各時間の0分のデータのみとする。

```py
import pandas as pd
import json
import csv
from pandas import json_normalize

#発電量に影響がある日中のデータ（8時から18時のみのデータ）とする
df_rain_analysis_q_short['date'] = df_rain_analysis_q_short["timestamp"].apply(lambda x : x.strftime('%Y/%m/%d'))
df_rain_analysis_q_short['time'] = df_rain_analysis_q_short["timestamp"].apply(lambda x : x.strftime('%H:%M:%S'))

df_rain_analysis_q_short_ = df_rain_analysis_q_short_\
    [(df_rain_analysis_q_short_['time']>="9:00:00") & (df_rain_analysis_q_short_['time']<="18:00:00")].reset_index(drop=True)

#0分のデータのみに抽出する
df_rain_analysis_q_short_['minites'] = df_rain_analysis_q_short_["timestamp"].apply(lambda x : x.strftime('%M:%S'))
df_rain_analysis_q_short_ = df_rain_analysis_q_short_\
    [df_rain_analysis_q_short_['minites'] == "00:00"]

```

出力結果は以下の通りです。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1194651/52e07ecf-fbb5-cb42-c1c6-c3fe850f2d35.png)


####予測と解析を結合する

~~~py
df_rain_analysis_q_short_ = df_rain_analysis_q_short_[["date", "rain_value"]]
df_rain_analysis_q_short_ = df_rain_analysis_q_short_.groupby('date').sum().reset_index()

df_forecast_ = df_forecast_0hour[['date', 'rain_forecast_0h']].copy()
df_forecast_ = df_forecast_.groupby('date').sum().reset_index()

total = pd.merge(df_forecast_, df_rain_analysis_q_short_, on='date', how = 'inner')
~~~

出力結果は以下の通りです。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1194651/c347c10a-bbb6-38db-1b9d-ccec73265cf2.png)


左のrain_dorecast_0hが0時時点での予測値、rain_valueが解析値となります。

2021/7/24では、予測では降雨量の予測は0に近かったにも関わらず、解析値が大きいことが分かります。
グラフで見ても、7/24は予測より非常に多い降雨量が発生していることが分かります。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1194651/d2f6facd-a698-b9c7-749d-4eeac1248cf8.png)


###2021/7/24 の太陽光発電量と解析雨量の関係を見てみる
####太陽光発電量と解析雨量のデータを結合する
・あらかじめ、太陽光発電量データを、「pv.csv」として配置する。

~~~py
import datetime

target_date = '2021-07-24'
target_date_next = '2021-07-25'

#太陽光発電量データの読み込み
df_pv = pd.read_csv("pv.csv")
df_pv['timestamp'] = pd.to_datetime(df_pv['timestamp'])
df_pv['timestamp'] = df_pv['timestamp'].dt.round("1H")
df_pv = df_pv.groupby('timestamp').sum().reset_index()

df_pv_ = df_pv[(df_pv['timestamp'] >= target_date) & (df_pv['timestamp'] < target_date_next)].reset_index(drop=True)

df_rain_analysis_q_short_['timestamp'] = \
    df_rain_analysis_q_short_['timestamp'].apply(lambda x : x.strftime('%Y/%m/%d %H:%M:%S'))
df_rain_analysis_q_short_['timestamp'] = pd.to_datetime(df_rain_analysis_q_short_['timestamp'])

#後ろラベルを前ラベルに変換
df_rain_analysis_q_short_['timestamp'] = df_rain_analysis_q_short_['timestamp'] + datetime.timedelta(hours=-1)
df_rain_analysis_q_short_ = df_rain_analysis_q_short_[["timestamp", "rain_value"]]

df_rain_analysis_q_short_ = df_rain_analysis_q_short_\
    [(df_rain_analysis_q_short_['timestamp'] >= target_date) & (df_rain_analysis_q_short_['timestamp'] < target_date_next)].reset_index(drop=True)

df_pv_ = pd.merge(df_pv_, df_rain_analysis_q_short_, on = "timestamp")
df_pv_ = pd.merge(df_pv_, df_forecast_0hour, on = "timestamp", how = "left")
~~~

####グラフで見てみる

~~~py
import numpy as np
import matplotlib.pyplot as plt

fig, ax1 = plt.subplots()

timestamp = df_pv_['timestamp']
pv = df_pv_['kwh']
rain_value = df_pv_['rain_value']
isolation_forecast_0h = df_pv_['isolation_forecast_0h']

ax2 = ax1.twinx()

ax1.plot(timestamp, pv, c="r", label="pv")
ax2.plot(timestamp, rain_value, c="b", label="rain_value")

ax1.set_xlabel("time")
ax1.set_ylabel("pv")
ax2.set_ylabel("rain_value")

labels = ax1.get_xticklabels()
plt.setp(labels, rotation=45, fontsize=10)

#凡例の表示
h1, l1 = ax1.get_legend_handles_labels()
h2, l2 = ax2.get_legend_handles_labels()
ax1.legend(h1 + h2, l1 + l2, loc='upper left', bbox_to_anchor=(1, 1), facecolor='gray')

ax1.grid()
~~~
9時に30mm以上の降雨量があり、その他はほぼ0mmであることから、最初に定義した内容でゲリラ豪雨であると言えます。
またゲリラ豪雨が発生した9時に、太陽光発電量が落ちていることが分かります。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1194651/f0abcf33-ebfa-4904-eba9-7b06a1a30c9f.png)

ゲリラ豪雨が発生した時間帯、太陽光発電量が瞬間的に落ち込むことが分かりました。
次回は、直前にゲリラ豪雨を検知して、太陽光発電量予測を修正することができるのか！？を
KCCS APIサービスの他のデータを使って、引き続き調査していきたいと思います。