Kugelblitz

いつ何時誰の挑戦でも受ける!

太陽高度と太陽方位角の算出

ちまちま進めているWebGLアプリの開発ですが、黄道(太陽の見かけ上の通り道)をリアルにしてみました。

デモページはここです。また、あまりスペックの良くないマシン用に低解像度版も作ってみました。

黄道をリアルに再現するためには、太陽高度と方位角を算出する必要がありますが、数式等は「■太陽の高度と方位角(JavaScriptによる簡易版)」を参考にしています。

詳細は、前述のページを見てもらえればだいたい分かるかと思いますが、かいつまんで説明すると、太陽高度と方位角を決める要素は、

  1. 太陽赤緯:太陽光線と地球の赤道面との角度
  2. 時角:天球上での天体の位置を表すのに用いられる値の一つ
  3. 経度:赤道面との角度

です。「時角」は馴染みのない言葉ですが、時間、って考えてもだいたいOKだと思います。

太陽高度と方位角を計算することで、日の出、日の入り、南中時間等が実際と同じになっています。

日時を決める

まず、何月何日の何時の太陽高度と方位角を算出したいのかを決めます。

var timer = +Date.now();
timer = timer / 100000.0;
// 元旦からの経過日数
var day = (timer % (24.0 * 365)) / 24.0; 
var hour = (timer % 24.0);

Date.now()で取ったエポック秒を100000で除算していますが、こうすることで、実際の1秒で、仮想空間上は1分の時間が経過するようにしています。リアル24分で仮想空間上1日ですね。この値を増やすと、仮想空間上の時間の進め方が速くなります。将来、ここをパラメータ指定できるようにするかもしれません。

あ、あとうるう年については考慮していませんね。

太陽赤緯・時角

元旦からの経過日数で、太陽赤緯を求めます。前述した参考サイトの数式をベタに実装しています。

        // 太陽赤緯
        var omega = Math.PI * 2 / 365; 
        var J = day + 0.5;
        var declination = (0.33281 - 22.984 * Math.cos(omega * J) - 0.34990 * Math.cos(2 * omega * J) - 0.13980 * Math.cos(3 * omega * J ) +
                3.7872 * Math.sin(omega * J) + 0.03250 * Math.sin(2 * omega * J) + 0.07187 * Math.sin(3 * omega * J)) *
                Math.PI / 180;

0.33281とかの値の意味はよくわかっていません^^。今度調べます。

次に、時角の計算です。

        // 時角
        var t1 = hour + (135 - 135) / 15;
        var hourAngle = (15 * t1 - 180) * Math.PI / 180;
        // 緯度
        var latitude = 35 * Math.PI / 180;        

135っていう数字がありますが、これは東経ですね。135は日本の東経です。(135 – 135) / 15 の結果は0なので、現状意味がないですが、将来的には東経の位置もパラメータ指定できるようにするために残してあります。

15っていう数字は、360度を24時間で割った数字で、1時間の差は、15度の経度の差になる、ってことです。

また、緯度は現状35度固定にしています。日本の緯度ですね。緯度を変えると白夜等を再現できるので、将来的にはここもパラメータ指定できるようにしたいです(逆に、経度はパラメータ指定できてもそんなに面白い結果にならないです)。

太陽高度・太陽方位角

材料はそろったので、いよいよ太陽高度・太陽方位角の算出です。

        //太陽高度
        var phi = Math.asin(Math.sin(latitude) * Math.sin(declination) +
                  Math.cos(latitude) * Math.cos(declination) * Math.cos(hourAngle));
        // 太陽方位角
        var sinA = Math.cos(declination) * Math.sin(hourAngle) / Math.cos(phi);
        var cosA = (Math.sin(phi) * Math.sin(latitude) - Math.sin(declination)) / Math.cos(phi) / Math.cos(latitude);
        var theta = Math.atan2(sinA, cosA) + Math.PI;

参考サイトの数式のまんまですが、一応数式の意味を理解した上で実装しました。理解していれば、sinとcosを間違えたとか、*と+を間違えたみたいなしょぼいバグをすぐに直せるようになります。

あと「均時差」っていう、天球上を一定な速さで動くと考えた平均太陽と、実際の太陽との移動の差があるらしいんですが、今回は実装していません。

球面座標系への変換

太陽高度・太陽方位角を求める手順は以上ですが、このアプリは太陽の位置を画面に表示してあげる必要があるので、太陽高度・太陽方位角を球面座標系に変換します。

        var radius = 450000;
        var x = (radius) * Math.cos(phi) * Math.cos(theta) ;
        var y = (radius) * Math.sin(phi);
        var z = (radius) * Math.cos(phi) * Math.sin(theta);

天球周りについては、次は雲の表現をちゃんとやってみたいと思います。

Pocket

他の記事