線分ＡＢとＸ軸が成す角θを取得[VB.NET][Math.Atan2][Parallel.For]
二点Ａ(x0, y0)、Ｂ(x1, y1)があったとする。
この時の線分ＡＢと水平に引かれたＸ軸が成す角θを求める。
Math.Atan2()関数を使うと簡単に求まります。
'点Aを原点として考える Dim px As Decimal = x1 - x0 Dim py As Decimal = y1 - y0 Dim ret As Double = Math.Atan2(py, px)
片方の点を原点として考える必要がある為、
点Ｂの座標から点Ａの座標分を引いて、
点Ａ、Ｂもろとも原点へずりっと移動してから、
Math.Atan2()に点Ｂの座標を与えます。
※注意！引数はY座標が第一引数です。お間違えなく。
今回はサンプル(Visual Studio2010 プロジェクト)を用意しました。
逆正接関数サンプルソース※レアです！今だけ！
以下が実行時のスクリーンショットです。
「他力本願」を選択して、画面上でマウスを３回（点Ａ決定、点Ｂ決定、確定）クリックすると、演算が始まるようになっています。
戻り値はラジアンで -π〜0〜πの範囲で返却されます。
画面ではそれを180 / πをかけて、度（°）に変換して表示しています。
これだけだとつまらないので、「Parallel.For」の評価も兼ねて、
逆正接の級数展開を利用して求めてみました。
逆正接関数の級数展開 Tan-1 x = x - 1/3 * x3 + 1/5 * x5 - 1/7 * x7 + 1/9 * x9 - ... Tanθ = py / px θ = Tan-1 (py / px) = (py / px) - 1/3 * (py / px)3 + 1/5 * (py / px)5 - 1/7 * (py / px)7 + 1/9 * (py / px)9 - ...
これを普通のfor文にすると、以下の様になります。
Dim p_nLoop As Int64 = 1000000
Dim theta As Decimal = 0D
Dim p As Decimal = py / px
For i As Int64 = 0 To p_nLoop - 1 Step 1
Dim k As Integer = (i + 1) * 2 - 1
Dim q As Decimal = 1 / k * Math.Pow(p, k)
If i Mod 2 = 0 Then
theta += q
Else
theta -= q
End If
Next
1000万回ループで2.2秒程度かかります。
※CPUは、Intel(R) Core(TM)2 Duo CPU E8500 @3.16GHz 3.17GHzです。
これを「Parallel.For」を使い、並列処理化させたのが次のfor文です。
Dim p_nLoop As Int64 = 1000000
Dim theta As Decimal = 0D
Dim p As Decimal = py / px
Parallel.For(0,
p_nLoop,
Function()
Return 0
End Function,
Function(i As Int64,
state As ParallelLoopState,
subtotal As Decimal)
If state.IsExceptional OrElse
state.IsStopped Then
'他でExceptionが発生した or Stopされた
state.Stop()
Return subtotal
End If
Dim k As Double = (i + 1) * 2 - 1
Dim q As Decimal = 0D
Try
q = 1 / k * Math.Pow(p, k)
Catch ex As Exception
'オーバーフローで計算できません！
Debug.WriteLine(ex.Message + " " + ex.StackTrace)
state.Stop()
Return subtotal
End Try
If i Mod 2 = 0 Then
subtotal += q
Else
subtotal -= q
End If
Return subtotal
End Function,
Sub(subtotal As Decimal)
theta += subtotal
End Sub)
1000万回ループで1.2秒程度になりました。
並列処理されていることの証拠に、タスクマネージャを見てみると、
CPU使用率の上限が50%だったのが、100%まで使われています:D
画面の「通常処理」が普通のfor文で、「並列処理」がParallel.Forでの処理です。
最終項は何回for文を回すか。最終項が多ければ多いほど、求める値は正確になります。
でも、今回わかったのは、素直にMath.Atan2()関数を使ったほうが速いということです。
「他力本願」での処理時間は、0.0000085秒前後です:-p
サンプルソースのプロジェクトはVisualStudio2010用ですが、
Parallel.For以外は特に2010特化の処理をしていないので、弄れば、2010以前でも動くと思います。
Parallel.ForはCPU100%張り付きにご注意下さい。
CPUが壊れた等の損害が生じても、責任は取れません。
参考文献
「新装版 オイラーの贈物 〜人類の至宝eiπ=-1を学ぶ〜」吉田武先生著
PLANEX MZK-WNH-BK
無線LANルータの調子が悪く、接続が途切れるようになったので、新調することにしました。
主に価格で以下を選びました。
PLANEX 150Mbps 高速無線LANルータ(簡易パッケージ) MZK-WNH-BK
http://www.amazon.co.jp/PLANEX-150Mbps-高速無線LANルータ-簡易パッケージ-MZK-WNH-BK/dp/B002WC9D1G/ref=sr_1_4?ie=UTF8&s=electronics
簡易パッケージという名の通り、化粧箱はなく、ちっこいダンボール箱に入ってました。
エコでいいと思います:-)
でもちゃんとマニュアルと設定ソフトCD付きなので、ご安心を。
機能的には前（IO-DATA WN-G54/R2）のものと差はなく、LANポートも4つ付いております。
接続までの解説…も必要ないです。
かんたんマニュアル通りに進めていけば、以下の設定完了画面に辿り着きます。
無線LANの設定は既に以下のような初期設定がされています。
SSID「planexuser」、WEPID「1」、WEPキーは本体底面に記載されているWEP128bitアスキー。
じゃあ、「機器の詳細設定」をしてみようということで、ボタンを押すと…。
出た。うおー、ログインIDとパスワードわかんないよーと混乱しないように。
デフォルトがID「admin」パスワード「password」って書いてあります…。
（ちなみに私は混乱してたので、てきとーにパスワードアタックして入ったのは内緒です）
設定画面は「http://192.168.1.1/」です。※デフォルト
最低限のものは揃っている感じです。うん、申し分ない。
PSPもDSも前回と同じ流れでちゃんと接続できました。うーん、最近のはかんたんだなぁ:-)
後は耐久性かな…。
WEBからHTMLを取得し、文字コードの自動判別までやってのけ
る。←あっ
WEBからHTMLを取得して、さらに文字判別までやってくれると便利ですよね。
自動判別には、DOBON.NET プログラミング道から「文字コードを判別する」を使用させて頂いてます。
神サイトです。
ただ、判別できない場合はNothingが返ってくるので、null判定が必要です。
Public Function GetHtml(ByVal url As String) As String
Try
Dim client As New Net.WebClient()
client.Headers.Add("User-Agent", "Mozilla/4.0 (compatible; MSIE 6.0; Windows XP)")
Dim srmBuff As IO.Stream = client.OpenRead(url)
Dim bytBuff As Byte() = {}
Using srmMemory As New IO.MemoryStream
Dim bytRead As Byte() = {}
Dim intRead As Integer = 0
Array.Resize(bytRead, 1024)
intRead = srmBuff.Read(bytRead, 0, bytRead.Length)
Do While intRead > 0
srmMemory.Write(bytRead, 0, intRead)
intRead = srmBuff.Read(bytRead, 0, bytRead.Length)
Loop
bytBuff = srmMemory.ToArray()
End Using
srmBuff.Close()
'文字コードを取得する
Dim enc As System.Text.Encoding = GetCode(bytBuff)
If enc Is Nothing Then
enc = System.Text.Encoding.ASCII
End If
'デコードして表示する
Dim html As String = enc.GetString(bytBuff)
Return html
Catch ex As Exception
Debug.WriteLine(ex.Message)
Throw ex
End Try
End Function
HTML解析をしたい場合に使えるかも知れません。
関数名から動的に関数呼び出しを行う
更新をさぼっているわけじゃないです。
ネタがないだけです。ネタがない＝障害がなく平和ってことなので、喜ばしい限りですが:-)
関数名が文字列として取得できていて、かつ、その関数を動的に呼び出したい
（コンパイル時にはどんな関数名かわからない）という場合、
次のようにして呼び出しを行うことができます。
'アセンブリ名WindowsApplication1のクラスMyClassに存在する関数method1を呼び出す。 'Imports System.Reflectionが必要です。 Dim classNm As String = "WindowsApplication1.MyClass" Dim methodNm As String = "method1" Dim asm As Assembly = Assembly.GetExecutingAssembly() Dim objClass As Object = asm.CreateInstance(classNm) Dim mi As MethodInfo = objClass.GetType.GetMethod(methodNm) mi.Invoke(objClass, Nothing)
簡単です。アセンブリ名は、プロジェクトのプロパティで決められるexe名です。
（判り易くするためにnull判定は省いていますので、各自お願いします。）
さらに引数付きの関数を呼び出す場合は以下のようにします。
Dim classNm As String = "WindowsApplication1.MyClass"
Dim methodNm As String = "method1"
Dim asm As Assembly = Assembly.GetExecutingAssembly()
Dim objClass As Object = asm.CreateInstance(classNm)
Dim mi As MethodInfo = objClass.GetType.GetMethod(methodNm)
Dim objList As Object() = {"aaa", "bbb", "ccc"}
mi.Invoke(objClass, objList)
MethodInfo::Invokeの第２引数にオブジェクト配列を渡します。
オブジェクト配列なので、構造体でもクラスでも引数として渡すことができます。
どんな時に使用するかといえば、
フレームワークを作ったりする場合に使うかもしれません。
設定ファイルから画面を構築したり…したくないですね。
クラスのプロパティにプロパティ名から動的にアクセスしたい場合は、
以下のようにします。
Imports System.Reflection
Public Class MyClass
Private _p1 As String
Public Sub New()
Me._p1 = String.Empty
End Sub
Public Property P1() As String
Get
Return Me._p1
End Get
Set(ByVal value As String)
Me._p1 = value
End Set
End Property
Public Sub SetValueWithPropertyName(ByVal propertyName As String, ByVal value As Object)
Dim obj As Object = Me
Dim pi As PropertyInfo = obj.GetType.GetProperty(propertyName)
If pi IsNot Nothing Then
Dim objValue As Object = value
pi.SetValue(obj, objValue, BindingFlags.Default, Nothing, Nothing, Nothing)
End If
End Sub
Public Function GetValueWithPropertyName(ByVal propertyName As String) As Object
Dim obj As Object = Me
Dim pi As PropertyInfo = obj.GetType.GetProperty(propertyName)
Dim objValue As Object = Nothing
If pi IsNot Nothing Then
objValue = pi.GetValue(obj, Nothing)
End If
Return objValue
End Function
End Class
小細工してクラスのメソッドにしちゃいましたが、
要はインスタンスからGetPropertyでPropertyInfoを取得して呼び出すだけです。
メソッドにした意味は特にありません。
Dim myClass As New MyClass() Dim value As String = String.Empty myClass.SetValueWithPropertyName("P1", "aaa") 'myClass.P1 = "aaa"と同じ value = myClass.GetValueWithPropertyName("P1") 'value = myClass.P1と同じ
やっぱりあまり使わないかもしれません。。。
平面直角座標変換（緯度経度→平面直角座標）(7)
３−９．ＸＹ座標を求める（その１）
いよいよ平面直角座標を求めます。
と、その前に必要な値を算出しておきます。
'# 経度の差（東方を正にとる） Dim deltaLambda As Decimal = radLng - radLng0 Dim eta2 As Decimal = e2 ^ 2 * Math.Cos(radLat) ^ 2 Dim t As Decimal = Math.Tan(radLat) '# 座標系の原点における縮尺係数 Dim m0 As Decimal = 0.9999 '# 卯酉線曲率半径 Dim w As Decimal = Math.Sqrt(1 - e1 ^ 2 * Math.Sin(radLat) ^ 2) Dim n As Decimal = a / w
３−１０．ＸＹ座標を求める（その２）
以上、算出された値から、ＸＹを求めます。
求められた値をテキストボックスに表示して終了します。
Dim x As Decimal = ((s - s0) + 1 / 2 * n * Math.Cos(radLat) ^ 2 * t * deltaLambda ^ 2 + _
1 / 24 * n * Math.Cos(radLat) ^ 4 * t * _
(5 - t ^ 2 + 9 * eta2 + 4 * eta2 ^ 2) * deltaLambda ^ 4 - _
1 / 720 * n * Math.Cos(radLat) ^ 6 * t * _
(-61 + 58 * t ^ 2 - t ^ 4 - 270 * eta2 + 330 * t ^ 2 * eta2) * _
deltaLambda ^ 6 - _
1 / 40320 * n * Math.Cos(radLat) ^ 8 * t * _
(-1385 + 3111 * t ^ 2 - 543 * t ^ 4 + t ^ 6) * deltaLambda ^ 8) * m0
Dim y As Decimal = (n * Math.Cos(radLat) * deltaLambda - _
1 / 6 * n * Math.Cos(radLat) ^ 3 * (-1 + t ^ 2 - eta2) * _
deltaLambda ^ 3 - _
1 / 120 * n * Math.Cos(radLat) ^ 5 * _
(-5 + 18 * t ^ 2 - t ^ 4 - 14 * eta2 + 58 * t ^ 2 * eta2) * _
deltaLambda ^ 5 - _
1 / 5040 * n * Math.Cos(radLat) ^ 7 * _
(-61 + 479 * t ^ 2 - 179 * t ^ 4 + t ^ 6) * deltaLambda ^ 7) * m0
txtX.Text = x.ToString("0.000#")
txtY.Text = y.ToString("0.000#")
ここまで読んでくださった方、お疲れ様でした。
以上で、緯度経度が平面直角座標に変換されるはずです。
平面直角座標変換（緯度経度→平面直角座標）(6)
３−６．赤道から指定緯度までの子午線弧長を求める（その１）
赤道から指定緯度までの子午線弧長を求めます。
求めるのは、赤道から原点緯度までの子午線弧長と、赤道から入力緯度までの子午線弧長の２つです。
まず、求める為のパラメータA〜Iまでを求めます。
Dim pA As Decimal = 1 + 3 / 4 * e1 ^ 2 + 45 / 64 * e1 ^ 4 + 175 / 256 * e1 ^ 6 + _ 11025 / 16384 * e1 ^ 8 + 43659 / 65536 * e1 ^ 10 + _ 693693 / 1048576 * e1 ^ 12 + 19324305 / 29360128 * e1 ^ 14 + _ 4927697775 / 7516192768 * e1 ^ 16 Dim pB As Decimal = 3 / 4 * e1 ^ 2 + 15 / 16 * e1 ^ 4 + 525 / 512 * e1 ^ 6 + _ 2205 / 2048 * e1 ^ 8 + 72765 / 65536 * e1 ^ 10 + _ 297297 / 262144 * e1 ^ 12 + 135270135 / 117440512 * e1 ^ 14 + _ 547521975 / 469762048 * e1 ^ 16 Dim pC As Decimal = 15 / 64 * e1 ^ 4 + 105 / 256 * e1 ^ 6 + 2205 / 4096 * e1 ^ 8 + _ 10395 / 16384 * e1 ^ 10 + 1486485 / 2097152 * e1 ^ 12 + _ 45090045 / 58720256 * e1 ^ 14 + 766530765 / 939524096 * e1 ^ 16 Dim pD As Decimal = 35 / 512 * e1 ^ 6 + 315 / 2048 * e1 ^ 8 + _ 31185 / 131072 * e1 ^ 10 + 165165 / 524288 * e1 ^ 12 + _ 45090045 / 117440512 * e1 ^ 14 + 209053845 / 469762048 * e1 ^ 16 Dim pE As Decimal = 315 / 16384 * e1 ^ 8 + 3465 / 65536 * e1 ^ 10 + _ 99099 / 1048576 * e1 ^ 12 + 4099095 / 29360128 * e1 ^ 14 + _ 348423075 / 1879048192 * e1 ^ 16 Dim pF As Decimal = 693 / 131072 * e1 ^ 10 + 9009 / 524288 * e1 ^ 12 + _ 4099095 / 117440512 * e1 ^ 14 + 26801775 / 469762048 * e1 ^ 16 Dim pG As Decimal = 3003 / 2097152 * e1 ^ 12 + 315315 / 58720256 * e1 ^ 14 + _ 11486475 / 939524096 * e1 ^ 16 Dim pH As Decimal = 45045 / 117440512 * e1 ^ 14 + 765765 / 469762048 * e1 ^ 16 Dim pI As Decimal = 765765 / 7516192768 * e1 ^ 16
…頭痛くなった方、同士です（笑
注意。VB.NETでは「^」は累乗を表す演算子ですが、C#.NETでは排他的論理和（xor）になってしまうので、
Math.Pow()関数を使うようにして下さい。
３−７．赤道から指定緯度までの子午線弧長を求める（その２）
求めたパラメータA〜Iよりさらにパラメータb1〜b9を求めます。
Dim bCoef As Decimal = a * (1 - e1 ^ 2) Dim b1 As Decimal = bCoef * pA Dim b2 As Decimal = bCoef * -pB / 2 Dim b3 As Decimal = bCoef * pC / 4 Dim b4 As Decimal = bCoef * -pD / 6 Dim b5 As Decimal = bCoef * pE / 8 Dim b6 As Decimal = bCoef * -pF / 10 Dim b7 As Decimal = bCoef * pG / 12 Dim b8 As Decimal = bCoef * -pH / 14 Dim b9 As Decimal = bCoef * pI / 16
３−８．赤道から指定緯度までの子午線弧長を求める（その３）
さて、やっとここでパラメータから子午線弧長を求めることができます。
お疲れ様でした^-^;
'# 赤道から座標系の原点の緯度までの子午線弧長
Dim s0 As Decimal = b1 * radLat0 + b2 * Math.Sin(2 * radLat0) + _
b3 * Math.Sin(4 * radLat0) + b4 * Math.Sin(6 * radLat0) + _
b5 * Math.Sin(8 * radLat0) + b6 * Math.Sin(10 * radLat0) + _
b7 * Math.Sin(12 * radLat0) + b8 * Math.Sin(14 * radLat0) + _
b9 * Math.Sin(16 * radLat0)
'# 赤道から入力緯度までの子午線弧長
Dim s As Decimal = b1 * radLat + b2 * Math.Sin(2 * radLat) + _
b3 * Math.Sin(4 * radLat) + b4 * Math.Sin(6 * radLat) + _
b5 * Math.Sin(8 * radLat) + b6 * Math.Sin(10 * radLat) + _
b7 * Math.Sin(12 * radLat) + b8 * Math.Sin(14 * radLat) + _
b9 * Math.Sin(16 * radLat)
平面直角座標変換（緯度経度→平面直角座標）(5)
３−３．ラジアンに変換
次に、先程得られた緯度経度１０進数を、ラジアンに変換します。
ラジアンというのは、180度で3.141592…（つまり、π）になるような数のことです。
半径１の円周上をどれだけ進んだかを表すことで、角度を表す方法なんですね。（弧度法といいます）
これは、度数にπ/180をかけると得られます。
180度でπとなる数なので、180 * π / 180 = π。うん、合ってる。
Dim radLat As Decimal = decLat * Math.PI / 180 Dim radLng As Decimal = decLng * Math.PI / 180 Dim radLat0 As Decimal = decLat0 * Math.PI / 180 Dim radLng0 As Decimal = decLng0 * Math.PI / 180
３−４．長半径、扁平率の決定
世界測地系か、日本測地系かで、使用する長半径、扁平率の値が異なります。
コンボボックスで選択された値によって、長半径、扁平率の値を変更します。
Dim a As Decimal = 6378137 '# 長半径（世界測地系） Dim f As Decimal = 1 / 298.257222101 '# 扁平率（世界測地系） If "日本測地系".Equals(cmbSokuchiKei.SelectedItem) Then a = 6377397.155 '# 長半径（日本測地系） f = 1 / 299.152813 '# 扁平率（日本測地系） End If
３−５．第一離心率、第二離心率を求める
ここらへんから、計算が激しくなります^-^;
離心率とは
意味はおいといて、計算式通りに計算していきます（笑
Dim e1 As Decimal = Math.Sqrt(2 * f - f ^ 2) '# 第一離心率 Dim e2 As Decimal = Math.Sqrt(2 * 1 / f - 1) / (1 / f - 1) '# 第二離心率
注意。VB.NETでは「^」は累乗を表す演算子ですが、C#.NETでは排他的論理和（xor）になってしまうので、
Math.Pow()関数を使うようにして下さい。
