海上の太陽光発電所に最適なとこってどこなんだろうか？

# はじめに
　最近ですが海上に太陽光発電所とか自然エネルギー関連の発電所を作ろうという動きが活発になっているようですね。最近テレビで見て、ふと「海の上に作るっていても、、どうやって場所を選定するんだろう」とちょっと気になりました。
　そこで何か調べてみたいなと思い、「太陽光発電といえば日射量だろう」と「海上に作るのだから水深も重要だろう」という私の偏見より、まずは日射量のヒートマップと水深のコンター図を重ねたマップを作って、「水深が比較的浅くかつ日射量が多い」という海上の太陽光発電所を建てるのに最適だと自分が思う場所を調査してみることにしました。

　本当は日射量だけじゃ絶対情報不足でしょうし、日射量が太陽光発電にどこまで影響するのかは私はかじったぐらいしか知りません。それから実際は波の高さとか船の航路じゃないかとか景観とか色々気にしなきゃいけないんでしょうけど、仕事じゃありませんし、そもそも趣味的に調べたいなと思ったことなのでそこはご了承ください。

# 実行環境
- CPU:  CORE i5 8th Gen
- メモリ: 8GB
- OSと環境: Window11にAnacondaをインストールした環境
- プログラミング言語とバージョン：python 3.11.5

# この記事で行ったこと
ステップ1. 日射量データのダウンロード(ERA5のダウンロード)
ステップ2. 日射量のヒートマップ作成
ステップ3. 水深データのダウンロード(ETOPO1のダウンロード)
ステップ4. 日射量のヒートマップと水深データのコンター図の重ね合わせ

## ステップ1. 日射量のダウンロード
- ダウンロードする気象データセット: ECMWF(欧州中期予報センター) ERA5
- ダウンロードする期間: ERA5で公開されている2022年の全ての日時
- ダウンロードする領域: 日本周辺　北緯:24°～46° 東経122°～153°（表示したい範囲は北緯:24°～45°だが余分に取得する）
- ダウンロードする項目: surface_thermal_radiation_downwards(地球表面の水平面に到達する太陽放射)
- ダウンロードする形式: netCDF

海上も含めた日射量のデータがどこにあるだろうと調べていたところ、参考サイト[1]のサイトよりECMWF(欧州中期予報センター)が公開しているERA5という気象データセットに行き当たりました。

　そしてERA5のデータのダウンロードの方法としては参考サイト[2]を参考にしました。 そのサイトの中で「自力で書いてもいいですが、ウェブサイトを利用すると便利です。ERA5の圧力面データの場合、下記のページで欲しいデータの情報を入力し、一番下の"Show API request"というタブをクリックすると、Pythonのプログラムが表示されます。」とのことでしたので、上記で紹介されている参考サイト[3]を利用して、以下のコードを作成して走らせデータをダウンロードしました。なお、その他の設定は上記のサイトをそのまま実行しました。

作成したダウンロード用コード
```
import cdsapi

c = cdsapi.Client()

c.retrieve(
    'reanalysis-era5-single-levels',
    {
        'product_type': 'reanalysis',
        'format': 'netcdf',
        'variable': 'surface_solar_radiation_downwards',
        'year': '2022',
        'month': [
            '01', '02', '03',
            '04', '05', '06',
            '07', '08', '09',
            '10', '11', '12',
        ],
        'day': [
            '01', '02', '03',
            '04', '05', '06',
            '07', '08', '09',
            '10', '11', '12',
            '13', '14', '15',
            '16', '17', '18',
            '19', '20', '21',
            '22', '23', '24',
            '25', '26', '27',
            '28', '29', '30',
            '31',
        ],
        'time': [
            '00:00', '01:00', '02:00',
            '03:00', '04:00', '05:00',
            '06:00', '07:00', '08:00',
            '09:00', '10:00', '11:00',
            '12:00', '13:00', '14:00',
            '15:00', '16:00', '17:00',
            '18:00', '19:00', '20:00',
            '21:00', '22:00', '23:00',
        ],
        'area': [
            46, 122, 24,
            153,
        ],
    },
    'download_surface_solar_radiation_downwards_2022.nc')
```
## ステップ2. 日射量のヒートマップ作成
年平均で出してもよかったのですが、月ごとのヒートマップを見たかったので「月ごとの平均」でヒートマップを作りました。以下は8月の日射量のヒートマップです。

![heatmap_month_8.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3581920/3c32d165-6de0-50b0-a97e-244c462a435f.png)


以下が作成したコードです。
```
import xarray as xr
import matplotlib.pyplot as plt
import cartopy.crs as ccrs

# NetCDFファイルの読み込み
data = xr.open_dataset('download_surface_solar_radiation_downwards_2022.nc')

# 日射量データの選択
solar_radiation = data['ssrd']
solar_radiation = solar_radiation/1000000

# 月ごとの平均を計算
monthly_mean = solar_radiation.groupby('time.month').mean('time')

# カラーバーの設定
colorbar_max = monthly_mean.max().values
colorbar_min = monthly_mean.min().values

# 各月の平均のヒートマップをプロット
for month in range(1, 13):
    plt.figure(figsize=(10, 8))  # アスペクト比を調整
    ax = plt.axes(projection=ccrs.PlateCarree())
    ax.set_extent([122, 153, 24, 45])  # 地図の表示範囲を設定
    ax.coastlines()  # 海岸線を描画

    # 月ごとのデータを選択
    monthly_data = monthly_mean.sel(month=month)

    # ヒートマップの作成
    plt.pcolormesh(data.longitude, data.latitude, monthly_data, shading='auto', transform=ccrs.PlateCarree(),vmin=colorbar_min, vmax=colorbar_max)
    plt.colorbar(label=f'Surface Solar Radiation Downwards in Month {month} (MJ/m²)')
    plt.title(f'Surface Solar Radiation Monthly Mean Heatmap in Month {month}')

    # PNGファイルとして保存
    plt.savefig(f"heatmap_month_{month}.png", bbox_inches='tight')
    plt.close()  # 現在の図を閉じる

```
## ステップ3. 水深データのダウンロード(ETOPO1のダウンロード)
- ダウンロードする地形データセット:ETOPO1
- ダウンロードしたデータファイル: ETOPO1_Bed_g_gmt4.grd

日射量データの次に必要な水深データのがどこにあるだろうと調べていたところ、参考サイト[4]よりETOPO1という地形データセットに行き当たりました。参考サイト[5]からブラウザ上でダウンロードできるのでダウンロード自体は簡単でした。解凍すると1GB近い大きなファイルで私の環境だと結構処理に時間がかかります。

## ステップ4. 日射量のヒートマップと水深データのコンター図の重ね合わせ
いよいよ、日射量のヒートマップと水深データのコンター図の重ね合わせですが、まずは結果がこちらです。結果は8月の日射量のヒートマップと水深のコンター図を重ね合わせました。水深は1000m毎に描いています。

![heatmap_waterdepth_month_8.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3581920/76a6df77-57ab-c964-7e7d-503837bb4d04.png)

コードについては、ETOPO1_Bed_g_gmt4.grdのデータが重いのでPCが落ちないように表示する範囲を絞ったり、普通に水深のプロットをしても左右上下が逆になったりと色々苦労しましたが、何とか修正したのが以下のコードです。
またステップ2では海岸線を描くのにax.coastlines() を使っていましたが、今回はETOPO1のデータを使って海岸線を描いたので滑らかになっています。

```
import xarray as xr
import matplotlib.pyplot as plt
import numpy as np
import cartopy.crs as ccrs
import rasterio
from rasterio.windows import from_bounds

# ETOPO1データの読み込み
with rasterio.open('ETOPO1_Bed_g_gmt4.grd') as etopo_dataset:
    window = from_bounds(122, 24, 153, 46, etopo_dataset.transform)
    etopo_data = etopo_dataset.read(1, window=window)
    etopo_transform = etopo_dataset.window_transform(window)
    etopo_x, etopo_y = np.meshgrid(
        np.linspace(etopo_transform.c, etopo_transform.c + etopo_transform.a * etopo_data.shape[1], etopo_data.shape[1]),
        np.linspace(etopo_transform.f, etopo_transform.f + etopo_transform.e * etopo_data.shape[0], etopo_data.shape[0])
    )

# NetCDFファイルの読み込み
data = xr.open_dataset('download_surface_solar_radiation_downwards_2022.nc')

# 日射量データの選択
solar_radiation = data['ssrd'] / 1000000 # MJ/m^2への変換

# 月ごとの平均の計算
monthly_mean = solar_radiation.groupby('time.month').mean('time')

# カラーバーの最大・最小値の定義
colorbar_max = monthly_mean.max().values
colorbar_min = monthly_mean.min().values

# 月ごとの平均のヒートマップのプロット
for month in range(1, 13):
    plt.figure(figsize=(12, 8))
    ax = plt.axes(projection=ccrs.PlateCarree())
    ax.set_extent([122, 153, 24, 45])

    # 月ごとの日射量データの選択
    monthly_data = monthly_mean.sel(month=month)

    # ヒートマップの作成
    mesh = ax.pcolormesh(data.longitude, data.latitude, monthly_data,
                         shading='auto', cmap='viridis',
                         vmin=colorbar_min, vmax=colorbar_max)

    # 水深のコンターを追加(1000m間隔)
    contour_levels = np.arange(-10000, 1, 1000)
    contours = ax.contour(etopo_x, etopo_y, etopo_data, levels=contour_levels,
                          colors='grey', linestyles='solid', linewidths=0.5, transform=ccrs.PlateCarree())
    coastline = ax.contour(etopo_x, etopo_y, etopo_data, levels=[0],
                           colors='black', linestyles='solid', linewidths=0.5, transform=ccrs.PlateCarree())

    # 水深値のラベル追加
    ax.clabel(contours, inline=True, fontsize=8, fmt='%1.0f m')

    # カラーバーの追加
    cbar = plt.colorbar(mesh, orientation='vertical', pad=0.02, aspect=50)
    cbar.set_label('Surface Solar Radiation Downwards (MJ/m²)', rotation=270, labelpad=20)

    # タイトルの設定
    plt.title(f'Surface Solar Radiation Monthly Mean Heatmap with Water Depth in Month {month}')

    # 図の保存
    plt.savefig(f"heatmap_waterdepth_month_{month}.png", bbox_inches='tight', dpi=300)
    plt.close()
```

# 考察と課題
これだけで考察するのは情報不足も甚だしいですが、やらせていただきます。

日射量が多くても電力消費地から何百キロも離れていれば、損失のほうが大きくなりそうです。それを考慮して考えると、適地として沖縄周辺、豊後水道や土佐湾よりもちょっと南あたりが日射量が大きく、かつ水深も比較的浅いように見えます。それから瀬戸内海も比較的日射量が大きいようにも見えます。

課題としては以下があげられます

- 調査期間が2022年の1年のみ
- 調査領域が日本周辺という大きい地域で大まかに適地を見ているだけ（より細かい解像度で見る必要あり）
- 調査対象が日射量と水深だけで、海岸からの距離、波高、潮の流れの速さ、台風の進路でないか、船の航路になってないか、など他の項目が不足している

これらすべてに取り組むのは難しいですが、何個かは今後取り組んでみたいと思います。

👆「調査領域が日本周辺という大きい地域で大まかに適地を見ているだけ（より細かい解像度で見る必要あり）」について[海上の太陽光発電所に最適なとこってどこなんだろうか？　候補地探索編](https://qiita.com/bebebefu/items/271fa3e2569e1dbea6d2)の記事で取り組んでみました。

👆[海上の太陽光発電所に最適なとこってどこなんだろうか？　波高調査編](https://qiita.com/bebebefu/items/19f887e7c8db56064551)の記事にて波高について調査してみました。

# 参考サイト
 [1] [ERA5の気象データを可視化する①～折れ線グラフ～](https://note.com/ats030/n/nfd8b40e1d066)
[2] [ERA5をダウンロードする ](https://biboroku-yk.blog.jp/archives/17262987.html)
[3] [The Climate Data Store:ダウンロード用のプログラムが表示されるサイト](https://cds.climate.copernicus.eu/cdsapp#!/dataset/reanalysis-era5-single-levels?tab=form)
[4] [ データはどこにあるのか？　海底地形図を描くためのデータを入手するには以下の４つの方法がある......](http://ofgs.aori.u-tokyo.ac.jp/~okino/gmtscripts/datainfo.html)
[5] [ETOPO Global Relief Modelのリンク](https://www.ncei.noaa.gov/products/etopo-global-relief-model)