2014年9月30日 (火)

Sunrise::データ表参照による計算

■データ表参照による計算
日出日入時刻等のデータの日々の変化は,毎年同じパターンの繰返しになる.
場所を固定すれば,基準年1年分の計算データから,補正により全ての年の
計算結果を導出できる.
事前計算データ
 基準年:2016年(通算日数の補正値が0になる)
 ファイルサイズ:2940バイト(4+367*8)
 データ内容:
   位置情報 2+2バイト
   日出時刻,日入時刻,日出方位,日入方位,南中高度 8バイト
 ファイル名:c:\_dat\srise$?$.zzd  (?に都市番号が入る)
 
処理内容
・良く利用する地点のデータを1年分保存する.
 データ保存機能を追加(10地点まで保存可能)
・利用の際は,都市番号と位置情報を照合し,一致した場合だけ使用する.
 不一致の場合,事前計算データがない場合は,通常の計算を行う.
・年ごとの補正は,整数計算で行い,浮動小数点は使わない.(高速)
・この計算は,一覧表画面の計算で用いる.
・メイン画面では,浮動小数点による近似式計算を用いる.
 
計算の比較
1.近似式による計算(==>部分がHP200LXでリアルタイムに行う計算)
  赤緯/均時差データ --> フーリエ係数計算 --> 赤緯/均時差計算 ==> 日出日入計算
2.データ表参照による計算
  赤緯/均時差データ --> 日出日入計算 ==> 年次変化補正計算
元になるデータが同じなので,両者の計算値には,ほとんど差がない.
差は全ての項目で,最大 0.1(端数処理の誤差程度).
 
[HP200LX/LXDIC]

|

2014年9月24日 (水)

Sunrise::近似式再考

■近似式再考
前回の結論は,通算日数に補正を加えれば,ある基準年の係数を用いて,
すべての年の赤緯/均時差の近似計算ができるということであった.
しかし,良く考えてみれば,わざわざ近似式を計算する必要がないことに
気付く.すなわち,赤緯/均時差の実データを直接利用し,補正を加えな
がら,他の年にも使い回わせば良い.
これなら,1件につき合計で 16回になる三角関数の計算が不要になり,
CPUの負荷が軽減される.特に,HP200LX のような非力なCPUを用いる場合
には好ましい.
 
ある年の赤緯/均時差の 366日分のデータをとり,これを元に補正を加え
ながら,他の全ての年に適用する.(次年の1/1を含むデータ)
(閏年を基準年とする場合の通算日数の補正)
 y%4==0 (閏年):  0    :J = days;
 y%4==1 の場合: +0.75 :J = days + 0.75;
 y%4==2 の場合: +0.5  :J = days + 0.5;
 y%4==3 の場合: +0.25 :J = days + 0.25;
J が整数ではなくなるので,求める値は,整数点のデータから,内挿によ
り計算する.
 
●データ量:2202バイト(6×367)(比較的小さなデータになる)
 赤緯 :long(4バイト)秒数で表し整数化
 均時差:int(2バイト) 秒数で表し整数化
(理科年表あるいは高精度計算サイトのデータを使うことが可能)
 
●計算結果
・4年分の赤緯/均時差データ表を用いた計算
・1年分の赤緯/均時差データ表(2016)を補正しながら使う計算
上記2つの計算の誤差:4秒程度(近似式を使うより誤差が小さい)
 2014: 日出 (-00.05, 00.02)分   日入 (-00.06, 00.03)分
 2015: 日出 (-00.03, 00.01)分   日入 (-00.04, 00.03)分
 2016: 日出 ( 00.00, 00.00)分   日入 ( 00.00, 00.00)分 <--基準年
 2017: 日出 (-00.03, 00.03)分   日入 (-00.02, 00.03)分
数年のスパンで考えれば,365.25日ごとに全く同一の変動パターンを繰り返す
と見なしても,誤差はほとんど生じない.
 
●三角関数の計算回数
 一覧表1日分: 67回 ==>  51回 ( -16回)
 グラフ1日分:1125回 ==> 1125回 (変らず)
 一覧表の都市切替時間 1.5秒 ==> 1.2秒 (7日分の計算時間)
 
●計算高速化の可能性
数年のスパンで考えれば,365.25日ごとに全く同一の変動パターンを繰り返す
のであるから,良く利用する都市の日出日入時刻等のデータについては,予め
1年分を計算しておけば,すべての年に対応できる.
浮動小数点計算の必要がなく,大幅な高速化が可能になる.
  
[HP200LX/LXDIC]

|

2014年9月22日 (月)

Sunrise 開発メモ

■Sunrise 開発メモ

●近似式の改良
赤緯/均時差の計算は近似式によっているが,これが全ての計算結果を左右
する.これまでの近似式はあまり精度が高くなく,日出日入時刻の計算誤差
が1~2分程度になっていた.
 
今回は,実データを与えて,より精度の高い係数を求めた.
この近似式は,赤緯/均時差を表す関数をフーリエ級数展開したもので,係
数は以下の計算で求められる.
 an = 1/π∫f(x)cos(nx)dx   (x=0 --> 2π)
 bn = 1/π∫f(x)sin(nx)dx   (x=0 --> 2π)
 f(x) = a0/2+Σ(an*cos(nx)+bn*sin(nx))  (n=1,2,...)
 
計算結果から,4次の係数までとれば,十分な精度が得られることが判明.
赤緯(計算結果)[単位:度]
 a0/2=  0.394472
  a1=-22.899438
  a2= -0.388787
  a3= -0.153999
  a4= -0.009749
  a5= -0.003751
  b1= +4.027593
  b2= +0.045014
  b3= +0.084584
  b4= +0.005171
  b5= +0.003806
均時差 (計算結果)[単位:時間]
 a0/2=  0.000334
  a1= +0.008998
  a2= -0.056444
  a3= -0.001506
  a4= -0.002352
  a5= -0.000155
  b1= -0.122402
  b2= -0.155564
  b3= -0.005427
  b4= -0.003044
  b5= -0.000344
 
●年変化の補正
赤緯/均時差の年変化は,実用上の精度に関しては,無いに等しい.
すなわち,365.25日を周期として同じ変動パターンを繰り返すと見るこ
とができる.
ただし,年ごとの変動パターンと暦の間には,毎年約 0.25日ずつのズレ
が生じるため,これを補正する必要がある.
 
1/1からの通算日数の計算に補正を加える.(yは西暦)
 y%4==0 (閏年): -0.75
 y%4==1 の場合:  0
 y%4==2 の場合: -0.25
 y%4==3 の場合: -0.50
毎年,変動パターンの同一位置より,約0.25日分前にずれて,1/1を基点
とする通算日数は大きくカウントされるため,その分を差引いて補正する.
閏年の翌年は,前年の閏日の挿入により,暦と変動パターンのずれがリ
セットされる.(正確には,0.03日分残る)
 
(例)閏年の1/1を基点にする場合の通算日数補正
 y%4==0 (閏年):  0    
 y%4==1 の場合: +0.75 
 y%4==2 の場合: +0.5  
 y%4==3 の場合: +0.25 
フーリエ係数は,基点となる年のデータで計算する.
 
(補足)
上記計算の誤差(1年365.2425日で計算)
 基準年:   0
 4年後:  +0.03(1.0-0.2425*4)
 8年後:  +0.06
 12年後:  +0.09
 16年後:  +0.12
 20年後:  +0.15

[HP200LX/LXDIC]

|

2014年9月17日 (水)

Sunrise 開発メモ

■Sunrise 開発メモ
追加メモ
・計算精度向上(時刻1分以内,角度1度以内)
・閏年対応(精度向上に伴い年変化が無視できない) 
・日付設定に"西暦"を追加
・”都市設定モード”を追加
 ノーマルフォントでファンクションガイドに都市名表示
 
年設定
Sr50lp1
閏年 
Sr50lp2
 
太陽方位の時間変化
Sr50b24h  
  
ノーマルフォントのファンクションガイド
Sr51bmn
  
[HP200LX/LXDIC]

|

2014年9月13日 (土)

Sunrise 開発メモ

■Sunrise 開発メモ
Ver 4.6c
 
現在の仕様 
 
●ウィンドウ構成
・メイン画面(1)
  数値表+グラフ1,又は,グラフ2
・設定ダイアログ(3)
  位置設定
  日付設定
  時刻設定
・データ一覧画面(1)
  7週分の日出日入等の数値表
・メニュー
  &Location, &Date, &Time, &Help, &Quit
  
  
●詳細内容
・グラフの種類
 グラフ1:太陽高度の時間変化+太陽方位の簡易表示
 グラフ2:太陽高度と方位の関係図
 (グラフ1)太陽高度の時間変化
  Sn46c02
 (グラフ2)太陽高度と方位の関係
  Sn46c01
 
・メイン画面のファンクションキー割当て
 F1  Help         'h'
 F2  再描画       '/'
 F3  Cities       '\'
 F4  夏時間       Del
 F5  現時刻       space
 F6  グラフ       right, left
 F7  一覧表       up,down
 F8   位置        'l'
 F9   時刻        't'
 F10  日付        'd'
 
 
・データ一覧画面のファンクションキー割当て
 F2  再描画        '/'
 F3  Cities        '\'
 F10 Cancel        ESC
(データ一覧画面)                           ※方位表示の不具合修正
 Sr46cerr
  
 
・登録できる都市数:10
 都市選択は,メイン画面と一覧表画面で数字キーにより行う.
 
・ヘルプ表示(F1キー)
・都市名一覧表示(F3キー)
 メイン画面への単なる上書表示で,簡易的なもの.
 消去するときは,F2で再描画
 (ヘルプ画面)表示中もすべての操作が有効
 Sn46c05
 (都市名表示)表示中もすべての操作が有効
 Sn46c03 
 
・データキャッシュによる高速化
 浮動小数点演算の負荷が大きく,計算に時間がかかるため,
 表示の遅延を回避する目的で以下のデータに実施する.
  日出日入一覧表用計算データ(50日分)
  グラフ用計算データ(10日分)
 
・設定ファイル
#(設定順序)
#   東経,北緯,時差,地名
#(設定項目)
#   東経:実数または整数  -180°<= lo <= 180°(20バイト以内)
#   北緯:実数または整数   -90°<= la <= 90°(20バイト以内)
#   時差:整数で指定  -12 <= h <= +12  (日本JST=UTC+9)
#         夏時間を考慮しない場合の時差を設定する.
#   地名:8 バイト以内.それ以上は無視される.
#(規則)
#   ・各項目は半角コンマ [,] で区切る.
#   ・各項目の間には,半角スペースを挿入して整形できる.
#   ・最大で10地点を設定できる.
#   ・11個目以降の設定は無視される.
#   ・#で始まる行は,コメント行
#   
139.7414, 35.6581,  +9, 東京
136.9167, 35.1667,  +9, 名古屋
135.0000, 35.0000,   9, 明石
127.6667, 26.2167,  +9, 那覇
141.3500, 43.0667,  +9, 札幌
-0.12806, 51.50778,  0, London
-73.54,   40.46,    -5, NewYork
 
・エラー表示
 設定ファイルエラーの表示(位置情報が全く取得できない場合)
 Sn46c06
 
[HP200LX/LXDIC]

|

2014年9月 7日 (日)

Sunrise Ver 3.8c

■Sunrise Ver 3.8c
開発メモ
 
追加
・緯度経度の制限を外す(全地域をカバー)
 東経: -180°~ +180°, 北緯: -90°~ +90°
・日出・日入時刻が存在しない場合が発生する.(高緯度の場合)
・日出日入で計算エラーが出ても,計算を中断せず高度と方位を計算
 1日中昼の場合は,日出00:00,日入24:00,昼長24:00 と表示
 1日中夜の場合は,日出--:--,日入--:--,昼長 0:00 と表示
・グラフ枠を3ドット下へ拡大
 高度0度の線が見えるようにする.
・グラフの現在時刻(又は指定時刻)の位置に縦線を入れる.
 
 
浮動小数点使用の影響(Borland C++ 3.1)
・三角関数を自由に使えるので,自在に計算できる.
・数学ライブラリをリンクするため,プログラムサイズが大きい.
・関数計算処理の負担が大きく,グラフ表示と一覧表示が遅くなる.
 
 
今後の課題
・場所選択方式の変更
 各都市を常時ファンクションキーに割当てると他の多くの機能の割当てができない.
 必要な時に選択可能な場所を表示し,選択する方式にする必要がある.
・計算速度の改善
 浮動小数点ではなく,整数による計算が可能かどうか検討する.
 計算が高速化すれば,純粋に表示に起因する部分の速度向上も必要になる.
 
 
実行例
高緯度の場合(北緯85度)
Sn38c1  
 
北極(北緯90度)
Sn38c6  
 
南半球(北緯 -35度)
Sn38c5  
 
[HP200LX/LXDIC]

|

2014年9月 5日 (金)

Sunrise Ver 3.7

■Sunrise Ver 3.7
 
変更点
・夏時間のON/OFFをできるようにする.
 ON :時差を+1時間にする.
 OFF:時差はそのまま.
・夏時間は一時的な設定で,保存されない.
・設定ファイルでの時差の設定は,夏時間を考慮しない.
 139.7414, 35.6581,  +9, 東京
 135.0000, 35.0000,  +9, 明石
 -0.12806, 51.50778,  0, London
 -73.54,   40.46,    -5, NewYork
・簡単なヘルプ表示機能を追加
 ヘルプファイル使用
 
表示例
通常の場合 時差-5時間
Snsmt01  
 
夏時間表示(UTC-4が反転表示)時差-4時間 [DELキーでON/OFF]
Snsmt02
 
ヘルプ表示(キー割当ては確定していない)
Snhelp
 
[HP200LX/LXDIC]

|

2014年9月 2日 (火)

太陽高度方位計算プログラム

■太陽高度方位計算プログラム
 Sunrise Ver 3.5
 
追加メモ 
・時差の設定を可能にする.
・夏時間への対応は,時差の設定により,手動で行う.
・設定ファイルの例
#   東経,北緯,時差,地名
139.7414, 35.6581,  +9, 東京
136.9167, 35.1667,  +9, 名古屋
135.0000, 35.0000,  +9, 明石
127.6667, 26.2167,  +9, 那覇
141.3500, 43.0667,  +9, 札幌
-0.12806, 51.50778, +1, London
-73.54,   40.46,    -4, NewYork
 
表示例
Sn_lndn2
 
Sn_tky
 
Sn_ny4
 
[HP200LX/LXDIC]

|

2014年8月31日 (日)

太陽高度方位計算プログラム

■太陽高度方位計算プログラム
 Sunrise Ver 3.2b
 
追加メモ 
・任意の緯度経度,日付,時刻の設定を可能にする.
・緯度経度を浮動少数点表示に変更
・テスト版では緯度経度の制限を小さくしている.日出日入が存在しない場合は,
 計算エラーがでる.(0除算,acos定義域エラー等)
・エラー処理により,計算を中断しハングアップを防ぐ.
 計算可能な場所&日付に切替えれば回復可能.
 
 
(表示例)
日付の設定画面 8/31 ==> 6/20 (最上部左側の1行が入力域)
Sun32b1
 
時刻設定画面 18:40 ==> 12:40
Sun32b2
 
経度緯度設定画面 139.7333==>135.0, 35.65==>25.502   
Sun32b3 
 
1項目の設定ごとに情報更新
Sun32b4
 
テスト版では経度に制限はない(ロンドンの計算値)
Sn32b24
 
●計算結果の比較
 
高精度計算サイト"keisan"
日の出日の入り(世界地名選択)
http://keisan.casio.jp/exec/system/1253955558
 
ロンドン -0.12806  51.50778(夏時間)海抜0m
--------------------------------------------------------------------------
西暦    日出   方位角  南中   高度   日入  方位角  時差
2014/08/31 06:11:12 74.9136 13:00:51 47.0532 19:49:25 284.7493 1
 
sunrise 3.2b の計算結果
日本時間  14:12   75   21:01   47   03:50   285
JST-8    6:12   75   13:01   47   19:50   285
 
概ね良好な結果
・時刻の誤差1~2分
・角度の誤差1度以内
  
(参考)
・タイムゾーンの名称:BST  イギリス夏時間
・協定世界時との時差:UTC+1
・日本時間との時差:JST-8
 
[HP200LX/LXDIC]

|

2014年8月27日 (水)

太陽高度計算プログラムの開発

■太陽高度計算プログラムの開発
 名称:Sunrize
 
●現在の状況 Sunrize Ver 2.7
・浮動小数点を用いた三角関数による計算
 誤差:日出日入:1~2分程度
    方位,高度:1度以内
・表示項目
  場所(緯度,経度)        
  日付
  日出時刻,日出方位
  南中時刻,南中高度
  日入時刻,日入方位
  太陽高度,太陽方位(現在時刻の値)
  グラフ表示(高度の時間変化,高度・方位関係)
  一覧表の表示
 
●実行画面例
 Sunr41
 
 Sunr42
 
 Sunr43
 
 
●特異な例(秋分点近傍の赤道上の値)
 太陽は,真東-真上-真西と動き,高度と時刻の関係はほぼ線型になる
 Sun031
  
 太陽の方位は,南中時刻の前後1時間で,ほぼ180度変化する
 Sun032
 
 赤道上では 8/27現在,太陽は北方向にある.(南中でなく北中)
 北(0度)を中心にグラフ上反時計回りに動く.西(280)<=北(0)<=東(90)
 Sun0add
 
 ※いずれも妥当な計算結果になっている.
 
[HP200LX/LXDIC]

|

«memo::太陽高度と方位の計算