日々適当

hibitekitou

メルカトル図法へ

xsi |2015-07-11

それじゃ、メルカトル図法に行ってみましょう。
メルカトル図法の式はシンプル。Wikipediaによると、地図中央の経度をλ0とし、与えられる経緯度がλ、φである場合、

x = λ - λ0

で、yは sinφ をtanhの逆関数にかけたもの、という事になるそうです。
tanhってなんぞや、となるわけですが、pythonなら大丈夫。mathモジュールを読み込めば、その理屈を知らずとも、tanhを適用すればいいだけで、逆関数もatanhとして用意されています(^^;

import math
app = Application

oObj = app.Selection(0)

oPntArr = oObj.ActivePrimitive.Geometry.Points.PositionArray

newPosArr = [[],[],[]]

lamda0 = 139.741389

for i in range(len(oPntArr[0])):
	x = oPntArr[0][i] - lamda0
	phy = math.radians( -oPntArr[2][i] )
	y = math.degrees( math.atanh( math.sin( math.fabs(phy) ) ) )
	if phy < 0:
		y = -y
	
	newPosArr[0].append(x)
	newPosArr[2].append(-y)
	newPosArr[1].append(0.0)

oObj.ActivePrimitive.Geometry.Points.PositionArray = newPosArr 


理屈は理解せずともw、簡単でいいですね。ただし、極点付近は極端に伸びるので、極点付近は削除の上、実行したほうがよいでしょう。



今回は80度線でぶった切ってみました。
そうすると、実行するとこんな風になります。



重ねてあるのはgoole mapから。低緯度の地域では一致しています。
高緯度の所は、なにせ



こんなふうにメッシュが引き伸ばされている為に、結果として高緯度側のメッシュ分割が粗く、精度が悪くなっております。

でも、対世界地図で活躍してくれそうですな。

コメント ( 0 )|Trackback ( )
 
コメント
 
コメントはありません。
コメントを投稿する
ブログ作成者から承認されるまでコメントは反映されません
 
名前
タイトル
URL
コメント
コメント利用規約に同意の上コメント投稿を行ってください。

数字4桁を入力し、投稿ボタンを押してください。