10 CLS : PRINT , , "Sun.bas -- SUNRISE and SUNSET": PRINT
11 ' $Id: sun.bas,v 1.4 2008/01/05 23:27:10 jno Exp $
20 PI = 3.14159: KEY OFF
30 RAD = 2 * PI / 360 'factor converts degrees to rads
40 DEG = 1 / RAD      'factor converts radians to degrees
50 DEF FNASIN (A) = ATN(A / ((1 - A * A) ^ .5))
60 DEF FNACOS (A) = -ATN(A / SQR(-A * A + 1)) + PI / 2
80 '.........................inputs.........................................
90 PRINT ,
100 INPUT ; "For location"; Locat$: PRINT : PRINT ,
101 INPUT ; "For axial inclination"; INC: PRINT : PRINT ,
102 INPUT ; "North latitude of location"; LAT: PRINT : PRINT ,
200 ' ..............dom.....month..............................
111 INPUT ; "For Month (Nnn)"; Month$: INPUT ; " DOM"; DOM: PRINT
601 IF Month$ = "Jun" THEN Days = ABS(DOM - 21)
602 IF Month$ = "Jul" THEN Days = DOM + 9
603 IF Month$ = "Aug" THEN Days = DOM + 9 + 31
604 IF Month$ = "Sep" THEN Days = DOM + 9 + 31 + 31
605 IF Month$ = "Oct" THEN Days = DOM + 9 + 31 + 31 + 30
606 IF Month$ = "Nov" THEN Days = DOM + 9 + 31 + 31 + 30 + 31
607 IF Month$ = "Dec" THEN Days = DOM + 9 + 31 + 31 + 30 + 31 + 30
608 IF Month$ = "Jan" THEN Days = 31 - DOM + 20 + 31 + 30 + 31 + 28
609 IF Month$ = "Feb" THEN Days = 28 - DOM + 20 + 31 + 30 + 31
610 IF Month$ = "Mar" THEN Days = 31 - DOM + 20 + 31 + 30
611 IF Month$ = "Apr" THEN Days = 30 - DOM + 20 + 31
612 IF Month$ = "May" THEN Days = 31 - DOM + 20
613 IF Days >= 182 THEN Days = 182 - ABS(Days - 182)
621 PRINT
120 'PRINT , : PRINT ; "    -- Solstice is on June 21--"
131 'PRINT ; "    Dec Jan Feb Mar Apr May Jun | Jun Jul Aug Sep Oct Nov Dec"
132 'PRINT ; "     9   31  28  31  30  31  20 -  9  31  31  30  31  30  21"
133 'PRINT ; "        171 140 112  81  51  20    9  40  71 101 132 162"
134 'PRINT
622 PRINT ; Month$; DOM; "(days after/before solstice:)"; Days
180 '........................declination.............................
190 DEC = INC * (COS(RAD * Days * 360 / 365.25))
220 '...............latitude - elevation of noon sun.............
240 IF LAT >= 90 OR LAT <= 0 THEN PRINT ; "LAT OUT OF LIMITS": GOTO 699
250 ELE = 90 - LAT + DEC
260 IF ELE <= 0 THEN PRINT , "Arctic Circle; no sun": GOTO 699
270 IF ELE > 90 THEN ELE = 180 - ELE
280 PRINT ; "Noon elevation of sun is"; ELE; "degrees."
350 '...................sunset /sunrise....................
380 LAT = RAD * LAT
381 DEC = RAD * DEC: SSAZICOS = SIN(DEC) / COS(LAT)
410 SSAZI = PI - FNACOS(SSAZICOS)
412 SUMAZI = DEG * SSAZI - 90
471 PRINT ; "sunrise at "; SUMAZI; "degrees No or So of East"
500 '...........times and day length .............
520 SINHR = SIN(SSAZI) / COS(DEC): TIME = FNASIN(SINHR) * DEG
540 IF DEC > 0 THEN TIME = 180 - TIME
550 DHR = TIME / 15: WHR = INT(DHR): EHR = 11 - WHR
560 M = INT(60 * (DHR - WHR)): LM = 60 - M
561 PRINT ; "Sunrise at approximately "; EHR; ":"; LM; "am,"
562 PRINT ; "Sunset  at approximately "; WHR; ":"; M; "pm LMT"
563 HOURS = 2 * WHR: MIN = 2 * M
564 IF MIN >= 60 THEN MIN = MIN - 60: HOURS = HOURS + 1
577 PRINT ; "Length of the day "; HOURS; "hours"; MIN; "minutes"
699 END