前回の続きです。
VSOLJ のメーリングリストでご指摘と解説を頂きました (vsolj 11105)。
ありがとうございます。
角度変換に asin とか atan は使っちゃダメで、
というかそれ以前の問題で、私アホなことに定義域とか値域のことも
すっかり忘れていました・・・orz
で、早速前回のコードをいじくりました。
今度はうまくいったように見えます。
次回は適当な星表を読み込んでみて、銀河座標でプロットしてみようと思います。
=======
RA <- 72
DEC<- -10
#-------------
aN <- 282.25
I <- 62.6
rad <- pi/180
deg <- 180/pi
RA <- RA*rad
DEC<- DEC*rad
aN <- aN*rad
I <- I*rad
#------------
x <- cos(DEC)*cos(RA-aN)
y <- sin(DEC)*sin(I) + cos(DEC)*sin(RA-aN)*cos(I)
z <- sin(DEC)*cos(I) - cos(DEC)*sin(RA-aN)*sin(I)
l <- atan2(y,x) *deg +33
l <- ifelse(l<0,360+l,l)
b <- atan2(z,sqrt(x*x+y*y)) *deg
cat("銀経=",sprintf("%3.2f",l),", 銀緯=",sprintf("%3.2f",b),"\n")
※コメント投稿者のブログIDはブログ作成者のみに通知されます