動作環境

Xeon E5-2620 v4 (8コア) x 2
32GB RAM
CentOS 6.8 (64bit)
openmpi-1.8.x86_64 とその-devel
mpich.x86_64 3.1-5.el6とその-devel
gcc version 4.4.7 (とgfortran)
NCAR Command Language Version 6.3.0
WRF v3.7.1を使用。
Python 3.6.0 on virtualenv

ベースコード geometry > 2つの緯度経度から距離を求める > Haversine formula使用

OsakaやLondonで一度ずれた時の距離を見てみた。

calc_distance_170331b.py

from math import sin, cos, radians, sqrt, asin

# codingrule:PEP8


def dist_on_sphere(pos0, pos1, radius=None):
    '''
    distance based on Haversine formula
    Ref: https://en.wikipedia.org/wiki/Haversine_formula
    '''
    if radius is None:
        radius = 6378.137  # km (Earth's radius)
    latang1, lngang1 = pos0
    latang2, lngang2 = pos1
    phi1, phi2 = radians(latang1), radians(latang2)
    lam1, lam2 = radians(lngang1), radians(lngang2)
    term1 = sin((phi2 - phi1) / 2.0) ** 2
    term2 = sin((lam2 - lam1) / 2.0) ** 2
    term2 = cos(phi1) * cos(phi2) * term2
    wrk = sqrt(term1 + term2)
    wrk = 2.0 * radius * asin(wrk)
    return wrk


def get_shifted_position(pos0, shift_lat_deg=None, shift_lng_deg=None):
    if shift_lat_deg is None:
        shift_lat_deg = 0.0
    if shift_lng_deg is None:
        shift_lng_deg = 0.0
    lat, lon = pos0
    return (lat + shift_lat_deg, lon + shift_lng_deg)

print(dist_on_sphere.__doc__)

# at Location1
Osaka1 = 34.702113, 135.494807
print("Osaka:")
Osaka2 = get_shifted_position(Osaka1, shift_lat_deg=1.0)
print("%.2f km for 1 deg latitude" % dist_on_sphere(Osaka1, Osaka2))
Osaka2 = get_shifted_position(Osaka1, shift_lng_deg=1.0)
print("%.2f km for 1 deg longitude" % dist_on_sphere(Osaka1, Osaka2))

# at Location2
print("London:")
London1 = 51.476853, 0.0
London2 = get_shifted_position(London1, shift_lat_deg=1.0)
print("%.2f km for 1 deg latitude" % dist_on_sphere(London1, London2))
London2 = get_shifted_position(London1, shift_lng_deg=1.0)
print("%.2f km for 1 deg longitude" % dist_on_sphere(London1, London2))

実行

$ python calc_distance_170331b.py 

    distance based on Haversine formula
    Ref: https://en.wikipedia.org/wiki/Haversine_formula
    
Osaka:
111.32 km for 1 deg latitude
91.52 km for 1 deg longitude
London:
111.32 km for 1 deg latitude
69.33 km for 1 deg longitude

v0.2

get_shifted_position()のデフォルト引数の処理が冗長であった。

calc_distance_170331b.py

from math import sin, cos, radians, sqrt, asin

# codingrule:PEP8


def dist_on_sphere(pos0, pos1, radius=None):
    '''
    distance based on Haversine formula
    Ref: https://en.wikipedia.org/wiki/Haversine_formula
    '''
    if radius is None:
        radius = 6378.137  # km (Earth's radius)
    latang1, lngang1 = pos0
    latang2, lngang2 = pos1
    phi1, phi2 = radians(latang1), radians(latang2)
    lam1, lam2 = radians(lngang1), radians(lngang2)
    term1 = sin((phi2 - phi1) / 2.0) ** 2
    term2 = sin((lam2 - lam1) / 2.0) ** 2
    term2 = cos(phi1) * cos(phi2) * term2
    wrk = sqrt(term1 + term2)
    wrk = 2.0 * radius * asin(wrk)
    return wrk


def get_shifted_position(pos0, shift_lat_deg=0.0, shift_lng_deg=0.0):
    lat, lon = pos0
    return (lat + shift_lat_deg, lon + shift_lng_deg)

print(dist_on_sphere.__doc__)

# at Location1
Osaka1 = 34.702113, 135.494807
print("Osaka:")
Osaka2 = get_shifted_position(Osaka1, shift_lat_deg=1.0)
print("%.2f km for 1 deg latitude" % dist_on_sphere(Osaka1, Osaka2))
Osaka2 = get_shifted_position(Osaka1, shift_lng_deg=1.0)
print("%.2f km for 1 deg longitude" % dist_on_sphere(Osaka1, Osaka2))

# at Location2
print("London:")
London1 = 51.476853, 0.0
London2 = get_shifted_position(London1, shift_lat_deg=1.0)
print("%.2f km for 1 deg latitude" % dist_on_sphere(London1, London2))
London2 = get_shifted_position(London1, shift_lng_deg=1.0)
print("%.2f km for 1 deg longitude" % dist_on_sphere(London1, London2))