2018年11月14日水曜日

緯度経度周りのUtil(プログラム)

最近利用頻度が高くなってきたので備忘録
1.緯度経度間の距離を求める
2.緯度経度の60進数表記から10進数表記へ変換

Excelで作業していたのでvbaです
評価検証用として作成(他より拝借)したので精度は検証できていませんがおおよそいけそうです
2018.12.15追記) 60 to 10変換にて各項目の桁数が少ない場合に計算が合わない件を修正
Excel上で利用する場合、この関数をモジュールに組み込んでから
=LatLng60To10(A1)などとします


Option Explicit

Const GRS80_A = 6378137#
Const GRS80_E2 = 6.69438002301188E-03
Const GRS80_MNUM = 6335439.32708317
 
Const PAI = 3.14159265358979
 
Function deg2rad(deg As Double) As Double
 
    deg2rad = deg * PAI / 180#
 
End Function

'国土地理院サイトにて計算結果を目視確認によりテスト
'https://vldb.gsi.go.jp/sokuchi/surveycalc/surveycalc/bl2stf.html
'***********************************************************************
' [説 明] 2地点間の直線距離を求める関数
' [引 数] lat1:地点1緯度 lng1:地点1経度 lat2:地点2緯度 lng2:地点2経度
' [戻り値] 距離差(m)
'***********************************************************************
Public Function calcHubenyToMeter(lat1 As Double, lng1 As Double, lat2 As Double, lng2 As Double) As Double
 
    Dim my As Double
    Dim dy As Double
    Dim dx As Double
    Dim sin As Double
    Dim w As Double
    Dim m As Double
    Dim n As Double
    Dim dym As Double
    Dim dxncos As Double
     
    my = deg2rad((lat1 + lat2) / 2#)
    dy = deg2rad(lat1 - lat2)
    dx = deg2rad(lng1 - lng2)
     
    sin = Math.sin(my)
    w = Math.Sqr(1# - GRS80_E2 * sin * sin)
    m = GRS80_MNUM / (w * w * w)
    n = GRS80_A / w
     
    dym = dy * m
    dxncos = dx * n * Math.Cos(my)
     
    'メートルで取得するので、キロメートルに変換し小数第一位までにする
    calcHubenyToMeter = Round(Math.Sqr(dym * dym + dxncos * dxncos), 2)
 
End Function

'LatLng60To10("35.39.30.959")
'35.6585997222222
'LatLng60To10("139.44.43.594")
'139.745442777778
'LatLng60To10("035.40.51.100")
'35.6808611111111

'LatLng60To10("35.39.09.959")
'35.6527663888889
'LatLng60To10("35.39.9.959")
'35.6527663888889

Public Function LatLng60To10(ByVal target As Variant) As Double

    Dim arr() As String
    Dim cal As Double
    Dim tmpSec As String
    Dim calsec As Double
    arr = Split(target, ".")

    If IsNull(target) Then
        Exit Function
    End If
    
    If UBound(arr) < 1 Then
        LatLng60To10 = ""
        Exit Function
    End If
    arr(0) = Right("00" & arr(0), 2)
    arr(1) = Right("00" & arr(1), 2)
    arr(2) = Right("00" & arr(2), 2)
    If UBound(arr) > 2 Then
        arr(3) = Right("000" & arr(3), 3)
        tmpSec = arr(2) & arr(3)
    Else
        tmpSec = arr(2) & "000"
    End If
    
    calsec = CDbl(tmpSec) / 1000
    
    cal = CDbl(arr(0)) + (CDbl(arr(1) / 60)) + CDbl(calsec / 3600)
    LatLng60To10 = CStr(cal)

End Function

0 件のコメント: