見出し画像

プログラミング系科目オール「C」の俺がPythonに一目惚れして局所天気図の描写に挑む話 ~第1話~

あらすじ

大学の必修科目に設定されたプログラミング(C言語)の講義を機に、一気にアンチプログラミング人間になってしまった俺ことうえはーす。とある講義で利用したPythonの使い勝手の良さに一目惚れし、AMeDAS観測点データから局所天気図を描き、熱的低気圧と雷雨性高気圧を表せないかと試みる。一応は描けたように見えた天気図であったが、どうやら本人はその出来に納得いかないようで……


環境について

話を進める前に、万が一にもこの記事を読んで似たようなことをしてみたいという人が現れた時のために、Pythonの開発環境について述べておく。

実行するためのプロンプトとプログラムを構築する開発環境には、それぞれ大学の講義の一環でインストールしたAnacondaとVisual Code Studioを利用している。ダウンロード方法や設定方法については詳しいサイトがゴロゴロあるはずなのでここでは割愛。

Pythonそのもののバージョン及び天気図の描写に用いたパッケージと使用したバージョンは以下の通り。一応、それぞれの用途についてもコメントアウトの形で記してしておく。

Python 3.8.18

Package             Version
------------------- ------------
numpy               1.24.4       # 数値計算
pandas              2.0.3        # CSVファイルの読み込み
scipy               1.10.1       # AMeDASや格子点の座標データ間を補間
matplotlib          3.7.2        # 図、グラフの描写
Cartopy             0.21.1       # 日本地図の描写

CSVファイルの確認

AMeDAS観測点の観測データについては以下の気象庁のサイトからダウンロードできる。今回は2018年7月15日の、気圧を観測している観測点に限定してダウンロード。突風の発生時間は15~16時だが追々時間変化も確認できるよう時間は09時から21時までとしておくか。

ダウンロードしたデータには様々な文字情報も付随している。できるだけPythonで読み込みやすくできるようにしてやりたい。

まず、CSVファイルを扱ったことのある人ならわかるだろうが、とにかく文字化け防止や読み込み時の形式指定が面倒くさすぎるのでデータは英数字だけにしてしまう。

また、本来ダウンロードしたCSVファイルにはデータを扱ううえで欠かせない情報の精度に関するデータ(品質情報)もある。しかしこのデータまで逐一読み込ませて参照していたのでは本来の予定を1年超過して、なおかつもう1年使うことになったノートPCくんの負担なので、泣く泣く消すことにする。

あと何より緯度経度のデータが無い。これなくして天気図など書きようがないので10進法表記で打ち込んでやる。そうしてできたCSVファイルはこんな感じ。

pref: 所在都道府県、site: 観測点名、lat: 緯度(latitude)、lon: 経度(longitude)、JSTXX: 日本時間XX時

所在地についても記録しているのは後々正確な位置にプロットできているかの確認のためだ。

日本地図の描写

天気図を描写するうえで日本地図は絶対に必要だ。まずはその描写からしていきたいが……ぶっちゃけてしまうと偉大なる先駆者の素晴らしいサイトに記載されていたコードをほとんどそのまま利用してしまっている。なので描写の仕方については当該サイトを参考されたし。

一応、描写範囲をあらかじめ指定しておくことで逐一書き換える必要の無いようにはした。少し先の話だが、格子点の作成にも使える。

# 描写範囲
W = 136 # 西端
E = 140 # 東端
S = 34 # 南端
N = 38 # 北端
-------------
# 一例としてこんな感じ
ax.set_extent([W, E, S, N]) # 図の範囲を指定するコード

# リスト[値, 値, …]にしておかないとエラー吐くぞ。何も知らなかった俺はいっぺんやらかした。

できた地図はこれだ。白地図では味気ないのでせめてもの海に色を付けてみた。

地図に観測点を表示する

まずはpandasでCSVファイルを読み込んでもらう。

import pandas as pd # pandasで読み込んでもらうためのおまじない

# CSV読み込み
df = pd.read_csv('gust_20180715.csv') # gust = 突風

このCSVファイルから緯度と経度のデータを抜き出したい。とりあえずこんな感じでやってみた。

# CSVから経度、緯度のデータを指定して抽出
lon = df['lon'] # 経度(x方向)。lonの列から抜き出してね、という意になるらしい。
lat = df['lat'] # 緯度(y方向)。latの(以下同文

print(lat)
print(lon) # 緯度から表示したほうが見やすいのでこの順にする。
---------------------------------------
0     34.628333
1     35.166667
2     35.400000
3     36.155000
4     34.940000
5     34.733333
6     34.068333
7     34.761667
8     35.275000
9     36.055000
10    36.588333
11    37.390000
12    36.708333
13    36.791667
14    37.106667
15    37.893333
16    36.661667
17    36.246667
18    36.045000
19    35.523333
20    34.753333
21    34.603333
22    34.975000
23    35.113333
24    34.603333
25    35.045000
26    36.405000
27    35.666667
28    35.438333
29    36.548333
30    35.990000
31    36.150000
32    35.691667
33    34.748333
34    34.123333
35    34.986667
36    37.488333
Name: lat, dtype: float64
0     137.093333
1     136.965000
2     136.761667
3     137.253333
4     136.580000
5     136.518333
6     136.193333
7     136.141667
8     136.243333
9     136.221667
10    136.633333
11    136.895000
12    137.201667
13    137.055000
14    138.246667
15    139.018333
16    138.191667
17    137.970000
18    138.108333
19    137.821667
20    137.711667
21    138.213333
22    138.403333
23    138.925000
24    138.841667
25    139.091667
26    139.060000
27    138.553333
28    139.651667
29    139.868333
30    139.073333
31    139.380000
32    139.750000
33    139.361667
34    139.520000
35    139.865000
36    139.910000
Name: lon, dtype: float64

とりあえずデータは正しく抜き出せているようだ。

ではここから日本地図にマーカーをプロットしてもらおう。

# AMeDAS観測点プロット
ax.scatter(lat, lon, marker="s", color="magenta", s=25, edgecolors="none", transform=ccrs.PlateCarree(), zorder=6)

# markerで形状指定。sはsquare、つまり四角。
# マーカーの形状によってはedgecolorsを指定するとエラー吐くことも。"+"はダメだった。
# transform=ccrs.PlateCarree()は正距円筒図法で描いてねというおまじない的なやつ。これ抜くと座標がずれたり表示されないことも。
# zorderは図の描写の優先順位。マーカープロットまでの段階であれば、とりあえずこいつの値を一番デカくしてやればよい。

plt.show()

よし、ここまで順調に来られたぞ。この調子でマーカーに地名も重ねて表示してやろう。……と意気込んでいたうえはーすであったが、実は初回作成時にここで大きく躓いている。次の話はこの辺りの話からしていこう。

【次の話】


この記事が気に入ったらサポートをしてみませんか?