日付と場所を入力して太陽の南中高度を求める JIS FULL BASIC 言語のプログラムを作ってみました。 公転周期を365日ちょうどにしていることや公転軌道を正円にしていることから、 若干の誤差があります。

110 INPUT PROMPT "何月？　":a
120 INPUT PROMPT "何日？　":b
130 INPUT PROMPT "北緯 （ 明石は３４.６度 ）？　":e
140 INPUT PROMPT "東経 （ 明石は１３５度 ）？　":f
160 IF a=1 THEN LET k=9
170 IF a=2 THEN LET k=9+31
180 IF a=3 THEN LET k=9+31+28
190 IF a=4 THEN LET k=9+31+28+31
200 IF a=5 THEN LET k=9+31+28+31+30
210 IF a=6 THEN LET k=9+31+28+31+30+31
220 IF a=7 THEN LET k=9+31+28+31+30+31+30
230 IF a=8 THEN LET k=9+31+28+31+30+31+30+31
240 IF a=9 THEN LET k=9+31+28+31+30+31+30+31+31
250 IF a=10 THEN LET k=9+31+28+31+30+31+30+31+31+30
260 IF a=11 THEN LET k=9+31+28+31+30+31+30+31+31+30+31
270 IF a=12 THEN LET k=9+31+28+31+30+31+30+31+31+30+31+30
280 IF a=12 THEN GOTO 510 ELSE GOTO 290
290 LET x=k+b
300 IF x<365/4 THEN LET y=90-(e+23.4)+23.4*(1-COS(2*PI*x/365)) ELSE GOTO 610
310 IF x<365/2 THEN LET h=720-146+146*4*x/365 ELSE GOTO 710
400 PRINT " "
410 PRINT "冬至からの日数： ";x;"日後"
420 PRINT "昼間の時間：　";INT(h*10/60+0.5)/10;"時間"
430 PRINT "太陽の南中高度： ";INT(y*10+0.5)/10;"度"
440 PRINT "太陽の南中時刻： 12時 ";INT(4*(135-f)*10+0.5)/10;"分"
450 GOTO 800
510 IF b<22 THEN GOTO 290 ELSE GOTO 520
520 LET x=b-22
530 GOTO 300
610 IF x<365/2 THEN LET y=90-e+23.4*SIN(2*PI*(x-365/4)/365) ELSE GOTO 630
620 GOTO 310
630 IF x<3*365/4 THEN LET y=90-(e-23.4)-23.4*(1-COS(2*PI*(x-365/2)/365)) ELSE GOTO 650
640 GOTO 310
650 LET y=90-e-23.4*SIN(2*PI*(x-3*365/4)/365)
660 GOTO 310
710 LET h=720+146-146*4*(x-365/2)/365
720 GOTO 400
800 END

　 CASIO COMPUTER CO., LTD. の高精度計算サイト （ http://keisan.casio.jp/ ） にいけば、 正確な太陽の南中高度がわかります。


　 ＜ Ｃ言語によるプログラム ＞

＜ JavaScript によるプログラム ＞

　何月 ？ ＝ 　月
　何日 ？ ＝ 　日
　北緯 ？ ＝ 　度　（ 明石は ３４．６ 度 ）
　東経 ？ ＝ 　度　（ 明石は １３５ 度 ）
　　
　冬至からの日数 ＝ 　日後
　昼間の時間　　　＝ 　時間
　太陽の南中高度 ＝ 　度
　太陽の南中時刻 ＝ 12 時 分

※ 参照： 大学生のための物理学 ＞ 天文学 ＞ 変動する太陽の南中時刻