科学と家事とプログラミング (python を中心に)

python 温度計測 湿度計測 DS18B20 USB9097

角度を扱うプログラムは、要注意なのです

天文計算や宇宙機の姿勢軌道制御では、角度の計算が欠かせませんが、 角度を扱う時は、色々と気が抜けません。 代表的には以下の3点。

  • 単位
  • 極性
  • 定義域 (周回の不確定性の扱いを含む)

単位

初等教育や日常生活では、deg を使う場面が多いけれど、プログラミング言語(C 言語や python)の 数値演算ライブラリの引数は rad に統一されています。 角度の単位の間違いを減らす一つの方法は、「deg 単位の変数は使わない」 という方針を禁欲的に守ることかな...と思っています。 例えば、引数で与えられた角度 [deg] と sin の値を表示するなら

import sys
import math
DEG = math.pi/180

x = float(sys.argv[1])*DEG # 変数は rad で保持
s = math.sin(x)
print("x %8.3f [deg]  sin(x) %6.3f [-]" % (x/DEG, s)) # 出力時に deg 単位に変換

x を deg 単位で保持した方が無駄が少なそう... といった、ケチな事を考えてはいけません。 ここでブレると、いずれ間違いが生じます。

もちろん、例外はあって、10deg 単位のメッシュを切りたいとか、 そういった場合には deg 単位の変数をその目的に限定して使うことはありますが、 角度一般の変数としては、必ず rad を使うことにしています。 プロジェクトや開発チームで統一されているとなお吉。

参考までに、print 出力の x/DEG という表記は、 x を deg 単位で表した値という意味を直接表現していて分かり良い表記と思っています。 入力側の x = 3*DEG も、1 deg が 3 つ分の角度ということで、こちらも 時々みかける deg2rad や rad2deg といった定数を乗ずる変換より、コードの 意味が明瞭と思います。このあたりは慣れや好みの問題かもしれませんけど。

定義域

定義域や周回の不確定性の扱いが問題なる例として、位相の区間の表現を取り上げます。 与えられた位相 x が位相区間 a, b に含まれるかどうかを判定することを考える。 単位円上の2点 a, b の位相角で位相の区間が与えられるものとする。 (問題を簡単にするため、ここでは a, b で定まる狭い側の円弧を採用する) 数直線であれば、a < b となるように選んで、a < x < b を判定すれば良いけれど、 相手が位相の場合には、そう簡単ではありません。 問題は -720, -360, 0, 360, -360, 720 [deg] といった位相角が、位相としては同じだということにつきます。 a, b の符号や大小関係を駆使して、場合分けによる判定式を構築できなくはないけど、 なかなか複雑になるので、出来れば避けたいところです。 ここでは、位相の定義域や a, b の数値的な大小関係に依存しない解法を考えてみます。 理由は、設計も検証も楽だから。

成す角で勝負

簡単に思いつくのは、a と b の中点 c をとって、c と x の成す角を計算する方法でしょうか。 c と a の成す角を基準値として判定できそうです。位相ではなくベクトルの成す角を 使って定式化しているので、位相の定義域や周回の不確定性には影響されません。 中点 c を計算するときは、かならずしも円周上で考える必用はなくて、ベクトル a, b の 単純平均(を規格化)で十分な場合が多いですね。

判定用に局所座標を導入する

他には、a または b を位相の原点に採った座標を導入して、その座標で判定することでしょうか。 よく設計されたベクトル行列演算ライブラリがあれば、こんな感じ

u = vector(az=a, el=0) # a 方向単位ベクトル
v = vector(az=b, el=0) # b 方向単位ベクトル

x = u.unit()       # a 方向を位相原点にとる
z = (u * v).unit() # * は外積
y = (z * x).unit() # y 方向 (位相=+π/2)
loc_fm_ref = mat(x, y, z)  # 基準座標 ==> 判定局所座標

v_loc = loc_fm_ref*v_loc
x_loc = loc_fm_ref*vector(az=x, el=0) # x 方向

# v_loc.az() と x_loc.az() の大小関係を調べる。具体的な判定規則は、
# .az() の仕様に依存する

計算量は多いけれど、a や b の定義域や極性をプログラム中で定義しなおす感じですね。 結果の良し悪しや限界も分かりやすいので、利用しやすい方法と思います。

定義域については、0~360deg と -180~180deg のどちらかという場合が多いけれど、 周回の不確定性を(部分的に)許容する場合もあるので、必ず確認が必用です。

極性

位相の測り方は、CCW が多いけれど、CW も無くはないので、明示されていないときは 必ず確認を。