Python」カテゴリーアーカイブ

Pythonによる勾配別標高グラフ完成

ちょっと前から作っていた勾配別の標高グラフを描画するプログラムができました。このグラフを作りたかった理由は、”道路勾配は常に一定ではない”ことにエントリで触れたかったからです。

坂道の場合はよく平均勾配という数値で厳しさを表現することがありますが、平均勾配の値が同じでも急勾配と緩勾配が交互に出てくるタイプのものもあれば、全体を通して大きな変化のないタイプもあります。

こういったタイプの違いは路線の個性とでもいうべき面白いポイントなので、客観的な数字で表現できる方法が欲しかったわけです。

開発過程で試験データとして使っていたのは篠ノ井にある茶臼山のデータです。1区間を100mに設定して表現するとこのような感じになります。

グラフ1本が距離100mに相当します。

平均すると6%前後という坂なのですが、中盤に10%近いセクションがいくつかあることや、後半は相対的に勾配が緩めであることが読み取れます。

一般的にほとんど平坦として認識されるルートについても描画をしてみました。例として採用したのは神奈川県の府中街道です。向ヶ丘遊園-久地付近までの区間を切り出しました。

ルートとしてはほとんど平坦ですが、途中に跨線橋が2ヶ所あります。そのため、実際に走ってみるとちょっとした坂が2つあるという認識になります。区間を50mに設定し、細かめに表現するとこんな感じです。

陸橋の部分に勾配が付いていることが読み取れます。こういった”実は坂がある”ルートの表現もしてみたかったので、満足のいく結果となりました。

制作の過程でPythonとBokehのことをまた学べたので、こちらはこちらで頭の整理もかねてエントリにしたいと考えています。

Pythonで勾配別標高グラフを作る

J Sportsのサイクルロードレース中継を見ていると坂の勾配が1km毎に表示された断面図が表示されることがあります。

ああいったものを作れないかと思い、Pythonとbokeh.plotの組み合わせで試行錯誤をしてみています。

現状ブログエントリなどで使っている断面図はGPSデータを散布図としてプロットし、線のように見せています。

勾配別の標高グラフはGPSデータを一定距離毎に積算して、前のポイントとの差異を取って勾配を出すことにしました。グラフの表現方法は棒グラフにしてみました。

サンプルとして使用しているのは茶臼山のデータです。500m毎でデータを区切り、区間の平均勾配を出しています。おおよそやりたいことは実現できたので、もう少し見た目を手直しして、コードをきれいにしたら改めてエントリにまとめてみたいと思っています。

今回試作したものです。高さで標高、色で勾配を表現しています。
従来型の断面図です。

Pythonで気象庁アメダスデータから風配図を作る

自転車に乗っていると気になるのは風の問題です。空気抵抗で走行に大きく影響が出るので、自転車で出かける日の風向・風速は非常に気になる情報です。

関東平野に対して長野県では風に悩まされることは少ないのではないかと思っていたのですが、実際に走ってみるとあまり変わらないか、状況次第では関東平野より風が強く感じることもあります。

しかしこれはあくまで感覚的な話なので視覚化して比較する方法を検討してみました。

風配図

特定の地点でどちらの方向からどれくらいの強さの風が吹いてくるかを表現する時に使うのが”風配図”です。パッと見はレーダーチャートに似ていますが、方位(出現頻度)と同時に風の強さの内訳を表現するのが特徴です。

方位・出現頻度・風速の内訳の3要素を1つの図で表現するので、ExcelやGoogleスプレッドシートの機能では難しそうに感じました。調べてみたところPythonには専用のライブラリがあるらしいことが分かったので、bokeh.plotの要領で作図できないか調べてみました。

Windroseライブラリ

風配図を描くのに使えるのがWindroseというライブラリです。

windrose’s documentation

[windrose]

インストールはPythonなのでpipを使います。

<python>
$ pip install windrose

他にmatplotlibとnumpyが必要になるので、必要に応じてこれらもインストールします。

インストールができたらサンプルのコードをコピー&ペーストして仕組みをチェックしてみます。基本的には風向(サンプルコード内ではwd)と風速(サンプルコード内ではws)を組にして読み込ませれば良いことが分かりました。

実際のコードと描画

サンプルコードで試した情報を踏まえて実際に気象庁のアメダスデータを読み込んで描画するコードを作成しました。元はサンプルの先頭にある風向出現率をパーセントで表すコードです。気象庁のアメダスデータに対応させるために風向の文字列を角度に読み替える処理と、欠測データをスキップする処理を追加しました。

スクリプトと同じディレクトリに”wind.csv”という名称で気象庁のアメダスデータ(風向・風速のみ)を配置してコードを実行すると風配図が描画されます。

<python>
import csv
from dataclasses import replace
from email import header
import string
from windrose import WindroseAxes
from matplotlib import pyplot as plt
import matplotlib.cm as cm
import numpy as np

#csvの読み取り
csv_file = open("./wind.csv", "r", encoding="ms932", errors="", newline="")
reader = csv.reader(csv_file)

#ヘッダ情報を5行分スキップ
i = 0
while i <= 5:
    header = next(reader)
    i = i + 1

#リストの宣言
ws = [] 
wd = []
wdfloat = []

#csvの値をリストに追加
for row in reader:
    #欠測データのある行はスキップ
    if len(row[1]) != 0 and len(row[3]) != 0:
        ws.append(float(row[1])) 
        wd.append(str(row[3])) 


#風向の文字列を角度の数値に変換
for strd in wd:
    dictitonary = {
        '北北東':'22.5',
        '東北東':'67.5', 
        '東南東':'112.5', 
        '南南東':'157.5', 
        '南南西':'202.5', 
        '西南西':'247.5', 
        '西北西':'292.5', 
        '北北西':'337.5', 
        '北東':'45', 
        '南東':'135', 
        '南西':'225', 
        '北西':'315', 
        '北':'0', 
        '東':'90', 
        '南':'180', 
        '西':'270',
        '静穏':'0'
        }

#風向文字列を辞書データで変換
    floatd = dictitonary[strd]
#float型にしてリストに追加
    wdfloat.append(float(floatd))

#風配図描画
ax = WindroseAxes.from_ax()
ax.bar(wdfloat, ws, normed=True, opening=0.8, edgecolor='white')
ax.set_legend()
plt.show()

今回は直近3ヶ月のアメダスデータ(時別)を使って描画を行いました。図にしてみると傾向は明瞭で、

  • 風向は北-北東と南西-西の出現率が高い
  • 特に西南西の風が吹く場合が多い
  • 風速は2.0m/s-5.9m/sが大半を占める

ことが読み取れます。やはり図にしてみると非常にわかりやすく感じます。

今回、データは期間の全データを使ってざっくり作図をしましたが、もっと期間を細かくしたり、あるいは時間帯別に集計してみるとまた違う傾向が見えてくると思います。何となく自分の中で”こうなのではないか?”と思っている仮説がいくつかあるので、それらの検証を行ってみたいと思っています。