TheBackShed.com

Forum Index : Microcontroller and PC projects : Sunrise-sunset times

Author

Message

TassyJim

Guru

Joined: 07/08/2011
Location: Australia
Posts: 6283

Posted: 04:54pm 11 Sep 2015

Copy link to clipboard

Print this post

I was asked to update an old MMBasic program from a few years ago.
Old enough to have line numbers.

This I have done and here it is.
I also went back and started from the original GWBasic (?) program, documenting the step as I went.

I will finish writing up the process as I think it will help others.
I can also improve the Variable report in MMEdit to make life easier.

There are lots of fantastic programs out there and they are crying to get rediscovered.

This offering works on Maximites running V4.5 as well as micromites running V4.7 and I see no reason why V4.6 would be a problem.

' Sunrise-Sunset
' TassyJim 12 Sep 2015
' Suitable for MMBasic V4.5 and above on Maximites and Micromites
'
' Latitude South and Longitude West are negative
' Time zones east of Greenwich (Australian) are negative
'
' Original program notes:
' This program by Roger W. Sinnott calculates the times of sunrise
' and sunset on any date, accurate to the minute within several
' centuries of the present. It correctly describes what happens in the
' arctic and antarctic regions, where the Sun may not rise or set on
' a given date. Enter north latitudes positive, west longitudes
' negative. For the time zone, enter the number of hours west of
' Greenwich (e.g., 5 for EST, 4 for EDT). The calculation is
' discussed in Sky & Telescope for August 1994, page 84.

IF MM.VER >4.05 THEN ' OPTION EXPLICIT introduced in V4.6
OPTION EXPLICIT
DIM Lat,Lon,Tz,D,M,Y
ENDIF

DO
INPUT "Lat, Long (deg)";Lat,Lon
INPUT "Time zone (hrs)";Tz
INPUT "Year, Month, Day";Y,M,D
IF Y<100 THEN Y = Y + 2000 ' allow for two digit years

PRINT suntimes$(Lat,Lon,Tz,Y,M,D)
PRINT ""
LOOP
END

FUNCTION suntimes$(Lat,Lon,Tz,Y,M,D)
' returns a string with sunrise and sunset times
LOCAL Ax(2), Dx(2), crlf$
LOCAL K1,Z0,G,D1,F,J,S,A,J3,T,TT,T0,n,L,V,U,W,P,B
LOCAL A5,A0,D0,DA,DD,A2,D2,H0,H1,H2,L0,L2,C0,C,H7
LOCAL Z,M8,W8,V0,V2,E,T3,H3,M3,N7,D7,AZ,D5,Z1,V1
LOCAL Sset$,Srise$,Msg1$,Msg2$,Msg3$,Msg4$

crlf$=CHR$(13)+CHR$(10)
K1=15*(PI/180)*1.0027379
Sset$="Sunset at "
Srise$="Sunrise at "
Msg1$="No sunrise this date"
Msg2$="No sunset this date"
Msg3$="Sun down all day"
Msg4$="Sun up all day"

Lon=Lon/360: Z0=Tz/24
G=1: IF Y<1583 THEN G=0
D1=INT(D): F=D-D1-.5
J=0-INT(7*(INT((M+9)/12)+Y)/4)
IF G<>0 THEN
S=SGN(M-9): A=ABS(M-9)
J3=INT(Y+S*INT(A/7))
J3=0-INT((INT(J3/100)+1)*3/4)
ENDIF
J=J+INT(275*M/9)+D1+G*J3
J=J+1721027+2*G+367*Y
IF F<0 THEN
F=F+1: J=J-1
ENDIF
T=(J-2451545)+F
TT=T/36525+1: ' TT = centuries from 1900
' LST at 0h zone time
T0=T/36525
S=24110.5+8640184.813*T0
S=S+86636.6*Z0+86400*Lon
S=S/86400: S=S-INT(S)
T0=S*360*PI/180
T=T+Z0
' Get Sun's Position
FOR n = 1 TO 2
' Fundamental arguments
' (Van Flandern &
' Pulkkinen, 1979)
L=.779072+.00273790931*T
G=.993126+.0027377785*T
L=L-INT(L): G=G-INT(G)
L=L*PI*2: G=G*PI*2
V=.39785*SIN(L)
V=V-.01000*SIN(L-G)
V=V+.00333*SIN(L+G)
V=V-.00021*TT*SIN(L)
U=1-.03349*COS(G)
U=U-.00014*COS(2*L)
U=U+.00008*COS(L)
W=-.00010-.04129*SIN(2*L)
W=W+.03211*SIN(G)
W=W+.00104*SIN(2*L-G)
W=W-.00035*SIN(2*L+G)
W=W-.00008*TT*SIN(G)
' Compute Sun's RA and Dec
S=W/SQR(U-V*V)
A5=L+ATN(S/SQR(1-S*S))
S=V/SQR(U):D5=ATN(S/SQR(1-S*S))
'R5=1.00021*SQR(U) ' we don't use R5 anywhere
Ax(n)=A5: Dx(n)=D5
T=T+1
NEXT n
T=T-1
IF Ax(2)<Ax(1) THEN Ax(2)=Ax(2)+PI*2
Z1=(PI/180)*90.833: ' Zenith dist.
S=SIN(Lat*(PI/180)): C=COS(Lat*(PI/180))
Z=COS(Z1): M8=0: W8=0: PRINT
A0=Ax(1): D0=Dx(1)
DA=Ax(2)-Ax(1): DD=Dx(2)-Dx(1)
FOR C0=0 TO 23
P=(C0+1)/24
A2=Ax(1)+P*DA: D2=Dx(1)+P*DD
' Test an hour for an event
L0=T0+C0*K1: L2=L0+K1
H0=L0-A0: H2=L2-A2
H1=(H2+H0)/2: ' Hour angle,
D1=(D2+D0)/2: ' declination,
' at half hour
IF C0<=0 THEN
V0=S*SIN(D0)+C*COS(D0)*COS(H0)-Z
ENDIF
V2=S*SIN(D2)+C*COS(D2)*COS(H2)-Z
IF SGN(V0)<>SGN(V2) THEN
V1=S*SIN(D1)+C*COS(D1)*COS(H1)-Z
A=2*V2-4*V1+2*V0: B=4*V1-3*V0-V2
D=B*B-4*A*V0
IF D >= 0 THEN
D=SQR(D)
IF V0<0 AND V2>0 THEN suntimes$=suntimes$ + Srise$
IF V0<0 AND V2>0 THEN M8=1
IF V0>0 AND V2<0 THEN suntimes$=suntimes$ + Sset$
IF V0>0 AND V2<0 THEN W8=1
E=(0-B+D)/(2*A)
IF E>1 OR E<0 THEN E=(0-B-D)/(2*A)
T3=C0+E+1/120: ' Round off
H3=INT(T3): M3=INT((T3-H3)*60)
suntimes$=suntimes$ + ftime$(H3,M3)
H7=H0+E*(H2-H0)
N7=-COS(D1)*SIN(H7)
D7=C*SIN(D1)-S*COS(D1)*COS(H7)
AZ=ATN(N7/D7)/(PI/180)
IF D7<0 THEN AZ=AZ+180
IF AZ<0 THEN AZ=AZ+360
IF AZ>360 THEN AZ=AZ-360
suntimes$=suntimes$ + ", azimuth "+STR$(AZ)+crlf$
ENDIF
ENDIF
A0=A2: D0=D2: V0=V2
NEXT
IF M8=0 AND W8=0 THEN ' Special msg?
IF V2<0 THEN suntimes$=suntimes$ + Msg3$+crlf$
IF V2>0 THEN suntimes$=suntimes$ + Msg4$+crlf$
ELSE
IF M8=0 THEN suntimes$=suntimes$ + Msg1$+crlf$
IF W8=0 THEN suntimes$=suntimes$ + Msg2$+crlf$
ENDIF
END FUNCTION

FUNCTION ftime$(h,m) ' a time format method that works on most versions of MMBasic
ftime$=RIGHT$(STR$(h+100),2)+":"+RIGHT$(STR$(m+100),2)
END FUNCTION

Jim

VK7JH
MMedit

isochronic
Guru

Joined: 21/01/2012
Location: Australia
Posts: 689

Posted: 06:12pm 11 Sep 2015

Copy link to clipboard

Print this post

Good stuff !

Are those last-of-sun's-disc on horizon times ?
The NOAA spreadsheet gives sunrise/sunset times but there is a disclaimer,
the values are not to be used for legalistic purposes. (Atmospheric refraction varies slightly with pressure etc) How about posting some comparison values ?

I guess MMBasic uses Microchips double precision trig functions ? I am wondering if say the equivalent code on Arduino etc would work nearly as well.

It would be good to see MMx being used for alternative / solar projects a lot more.

TassyJim

Guru

Joined: 07/08/2011
Location: Australia
Posts: 6283

Posted: 08:48pm 11 Sep 2015

Copy link to clipboard

Print this post

The program uses

"Z1=(PI/180)*90.833: ' Zenith dist. "

So I think that 90.833 means last-of-sun's-disc horizon times.

As far as accuracy, MMBasic uses 32bit floating point and so did GWBasic/Quick Basic.

The times seem to agree with most on-line and BOM figures.
I deliberately didn't use seconds in the output.

Most of the variations between civil, military etc are for twilight time.

Jim

VK7JH
MMedit

MMAndy
Regular Member

Joined: 16/07/2015
Location: United States
Posts: 91

Posted: 05:28pm 12 Sep 2015

Copy link to clipboard

Print this post

Using your "old" sunrise sunset algorithm, in DOS, I inputted the following for Columbus, Ohio USA:

Lat: 39.9833 <--- in decimal
Long: -83.9833
timezone: 4 (EDT)

Year 2015
Day 12
Month 9

Old result: +- 1 min <--- good
sunrise at 07:11, azimuth 84.1997 <--- OK
sunset at 19:45, azimuth 275.545 <--- OK

New algorithm result for MM+:
sunrise at 03:02, azimuth 118.191 <---?
sunset at 12:35, azimuth 241.728 <---?

Did you change the input format?

TassyJim

Guru

Joined: 07/08/2011
Location: Australia
Posts: 6283

Posted: 06:05pm 12 Sep 2015

Copy link to clipboard

Print this post

Using your data:

DOS MMBasic Version 4.5
Copyright 2011-2014 Geoff Graham

Lat, Long (deg)? 39.9833,-83.9833
Time zone (hrs)? 4
Year, Month, Day? 15,9,12

Sunrise at 07:14, azimuth 83.8052
Sunset at 19:50, azimuth 275.939

I think you had the date order wrong.
I used Year, Month, Day to get away from the usual problem of some countries have their date the wrong way around. ie MM/DD/YY instead of the Australian DD/MM/YY

If you program gets it's input from a user, you will need to put sanity checks in to be confident that the data inputted is OK.
Even then, it is difficult to know if 1/4 is January 4th or 1st April.

The function doesn't do ANY range checking.

It is always possible that there is an error in my conversion.
I did the checking using my location which is South and West rather than your North and East.

Jim
Edited by TassyJim 2015-09-14

VK7JH
MMedit

isochronic
Guru

Joined: 21/01/2012
Location: Australia
Posts: 689

Posted: 07:22pm 12 Sep 2015

Copy link to clipboard

Print this post

The reason I ask, is that [I think] the trig libraries called in C default to industrial-strength high precision algorithms, even if the calling input and output is by using single precision variables. As an exercise I changed a program using dp variables to sp, with very little effect on the results ! The downside is the libraries are larger and a little slower, less precision routines can be specified if wanted but it is a deliberate step. Arduino "double" actually is single, so I am wondering what libraries arduino uses...may be it has space restraints more. There seem to be quite a few solar projects with arduino but not much detail on how useful the results are.

edit BTW why aren't the "number of posts" numbers increasing ? Perhaps the corgis are getting sneaking bytes

Edited by chronic 2015-09-14

TassyJim

Guru

Joined: 07/08/2011
Location: Australia
Posts: 6283

Posted: 07:46pm 12 Sep 2015

Copy link to clipboard

Print this post

I 'think' that MMBasic uses single precision internals for the trig functions.

I was raised with a slide-rule and lookup tables. Double precision wasn't an issue.

I do remember having to calculate the effects that small input errors would have on the end result. It is rather frightening sometimes.
I also should look at the processing order to see if there is any unnecessary loss of precision. That's more an issue with integer maths.

As far as this program goes, most of the inaccuracy will be from the ' Fundamental arguments - (Van Flandern & Pulkkinen, 1979)
Some of them are 9 significant digits long even though the floating point number they are stored in will round them off. Other constants are only 3 or 4 significant digits.

Within a minute is good enough for me and most solar tracking purposes.
I live on high ground so I get a head start....

The main thing is to know what precision your purposes require.

Jim

VK7JH
MMedit

TassyJim

Guru

Joined: 07/08/2011
Location: Australia
Posts: 6283

Posted: 09:58pm 12 Sep 2015

Copy link to clipboard

Print this post

I found a bug...
It's one that always seems to catch me out. Variables passed to functions are passed by reference so if your function changes it, the 'master copy' is also changed.
If you were to run my function twice in a row without resetting the variables, you would get different results.
The suntimes$() function now doesn't alter any of the variables passed to it.

This version also has a function to print out the inputs in a different format so the user can see that his/her input is correct (or not).

' Sunrise-Sunset
' TassyJim 13 Sep 2015
' Suitable for MMBasic V4.5 and above on Maximites and Micromites
'
' Latitude South and Longitude West are negative
' Time zones east of Greenwich (Australian) are negative
'
' Original program notes:
' This program by Roger W. Sinnott calculates the times of sunrise
' and sunset on any date, accurate to the minute within several
' centuries of the present. It correctly describes what happens in the
' arctic and antarctic regions, where the Sun may not rise or set on
' a given date. Enter north latitudes positive, west longitudes
' negative. For the time zone, enter the number of hours west of
' Greenwich (e.g., 5 for EST, 4 for EDT). The calculation is
' discussed in Sky & Telescope for August 1994, page 84.

IF MM.VER >4.05 THEN ' OPTION EXPLICIT introduced in V4.6
OPTION EXPLICIT
DIM Lat,Lon,Tz,D,M,Y
ENDIF

DO
INPUT "Lat, Long (deg)";Lat,Lon
INPUT "Time zone (hrs)";Tz
INPUT "Year, Month, Day";Y,M,D
IF Y<100 THEN Y = Y + 2000 ' allow for two digit years
PRINT ""
PRINT checkinput$(Lat,Lon,Tz,Y,M,D)
PRINT ""
PRINT suntimes$(Lat,Lon,Tz,Y,M,D)
PRINT ""

LOOP
END

FUNCTION suntimes$(Lat,Lon,Tz,Y,M,D)
' returns a string with sunrise and sunset times
LOCAL Ax(2), Dx(2), crlf$, Lonx, Dy
LOCAL K1,Z0,G,D1,F,J,S,A,J3,T,TT,T0,n,L,V,U,W,P,B
LOCAL A5,A0,D0,DA,DD,A2,D2,H0,H1,H2,L0,L2,C0,C,H7
LOCAL Z,M8,W8,V0,V2,E,T3,H3,M3,N7,D7,AZ,D5,Z1,V1
LOCAL Sset$,Srise$,Msg1$,Msg2$,Msg3$,Msg4$
Dy=D
crlf$=CHR$(13)+CHR$(10)
K1=15*(PI/180)*1.0027379
Sset$="Sunset at "
Srise$="Sunrise at "
Msg1$="No sunrise this date"
Msg2$="No sunset this date"
Msg3$="Sun down all day"
Msg4$="Sun up all day"

Lonx=Lon/360: Z0=Tz/24
G=1: IF Y<1583 THEN G=0
D1=INT(Dy): F=Dy-D1-.5
J=0-INT(7*(INT((M+9)/12)+Y)/4)
IF G<>0 THEN
S=SGN(M-9): A=ABS(M-9)
J3=INT(Y+S*INT(A/7))
J3=0-INT((INT(J3/100)+1)*3/4)
ENDIF
J=J+INT(275*M/9)+D1+G*J3
J=J+1721027+2*G+367*Y
IF F<0 THEN
F=F+1: J=J-1
ENDIF
T=(J-2451545)+F
TT=T/36525+1: ' TT = centuries from 1900
' LST at 0h zone time
T0=T/36525
S=24110.5+8640184.813*T0
S=S+86636.6*Z0+86400*Lonx
S=S/86400: S=S-INT(S)
T0=S*360*PI/180
T=T+Z0
' Get Sun's Position
FOR n = 1 TO 2
' Fundamental arguments
' (Van Flandern &
' Pulkkinen, 1979)
L=.779072+.00273790931*T
G=.993126+.0027377785*T
L=L-INT(L): G=G-INT(G)
L=L*PI*2: G=G*PI*2
V=.39785*SIN(L)
V=V-.01000*SIN(L-G)
V=V+.00333*SIN(L+G)
V=V-.00021*TT*SIN(L)
U=1-.03349*COS(G)
U=U-.00014*COS(2*L)
U=U+.00008*COS(L)
W=-.00010-.04129*SIN(2*L)
W=W+.03211*SIN(G)
W=W+.00104*SIN(2*L-G)
W=W-.00035*SIN(2*L+G)
W=W-.00008*TT*SIN(G)
' Compute Sun's RA and Dec
S=W/SQR(U-V*V)
A5=L+ATN(S/SQR(1-S*S))
S=V/SQR(U):D5=ATN(S/SQR(1-S*S))
'R5=1.00021*SQR(U) ' we don't use R5 anywhere
Ax(n)=A5: Dx(n)=D5
T=T+1
NEXT n
T=T-1
IF Ax(2)<Ax(1) THEN Ax(2)=Ax(2)+PI*2
Z1=(PI/180)*90.833: ' Zenith dist.
S=SIN(Lat*(PI/180)): C=COS(Lat*(PI/180))
Z=COS(Z1): M8=0: W8=0: PRINT
A0=Ax(1): D0=Dx(1)
DA=Ax(2)-Ax(1): DD=Dx(2)-Dx(1)
FOR C0=0 TO 23
P=(C0+1)/24
A2=Ax(1)+P*DA: D2=Dx(1)+P*DD
' Test an hour for an event
L0=T0+C0*K1: L2=L0+K1
H0=L0-A0: H2=L2-A2
H1=(H2+H0)/2: ' Hour angle,
D1=(D2+D0)/2: ' declination,
' at half hour
IF C0<=0 THEN
V0=S*SIN(D0)+C*COS(D0)*COS(H0)-Z
ENDIF
V2=S*SIN(D2)+C*COS(D2)*COS(H2)-Z
IF SGN(V0)<>SGN(V2) THEN
V1=S*SIN(D1)+C*COS(D1)*COS(H1)-Z
A=2*V2-4*V1+2*V0: B=4*V1-3*V0-V2
Dy=B*B-4*A*V0
IF Dy >= 0 THEN
Dy=SQR(Dy)
IF V0<0 AND V2>0 THEN suntimes$=suntimes$ + Srise$
IF V0<0 AND V2>0 THEN M8=1
IF V0>0 AND V2<0 THEN suntimes$=suntimes$ + Sset$
IF V0>0 AND V2<0 THEN W8=1
E=(0-B+Dy)/(2*A)
IF E>1 OR E<0 THEN E=(0-B-Dy)/(2*A)
T3=C0+E+1/120: ' Round off
H3=INT(T3): M3=INT((T3-H3)*60)
suntimes$=suntimes$ + ftime$(H3,M3)
H7=H0+E*(H2-H0)
N7=-COS(D1)*SIN(H7)
D7=C*SIN(D1)-S*COS(D1)*COS(H7)
AZ=ATN(N7/D7)/(PI/180)
IF D7<0 THEN AZ=AZ+180
IF AZ<0 THEN AZ=AZ+360
IF AZ>360 THEN AZ=AZ-360
suntimes$=suntimes$ + ", azimuth "+STR$(AZ)+crlf$
ENDIF
ENDIF
A0=A2: D0=D2: V0=V2
NEXT
IF M8=0 AND W8=0 THEN ' Special msg?
IF V2<0 THEN suntimes$=suntimes$ + Msg3$+crlf$
IF V2>0 THEN suntimes$=suntimes$ + Msg4$+crlf$
ELSE
IF M8=0 THEN suntimes$=suntimes$ + Msg1$+crlf$
IF W8=0 THEN suntimes$=suntimes$ + Msg2$+crlf$
ENDIF
END FUNCTION

FUNCTION ftime$(h,m) ' a time format method that works on most versions of MMBasic
ftime$=RIGHT$(STR$(h+100),2)+":"+RIGHT$(STR$(m+100),2)
END FUNCTION

FUNCTION checkinput$(Lat,Lon,Tz,Y,M,D)
LOCAL mydate$, mylat$, mylon$, NS$, de, mm, La, Lo
IF M>0 AND M<=12) THEN
mydate$=STR$(D)+"-"+MID$("xxJanFebMarAprMayJunJulAugSepOctNovDec",M*3,3)+"-"+STR$(Y)
ELSE
mydate$="Unknown date format"
ENDIF
IF D < 1 OR D >31 THEN
mydate$="Unknown date format"
ENDIF
IF Lat<0 THEN
NS$="South"
ELSE
NS$="North"
ENDIF
La = ABS(Lat)
de=INT(La): mm = INT((La - de)*60 + 0.5)
mylat$=STR$(de)+"deg "+STR$(mm)+"min "+NS$
IF Lon<0 THEN
NS$="East"
ELSE
NS$="West"
ENDIF
Lo = ABS(Lon)
de=INT(Lo): mm = INT((Lo - de)*60 + 0.5)
mylon$=STR$(de)+"deg "+STR$(mm)+"min "+NS$
IF Tz < 0 THEN
NS$="Ahead of"
ELSE
NS$="Behind"
ENDIF
checkinput$=mylat$+", "+mylon$+CHR$(13)+CHR$(10)
checkinput$=checkinput$+"("+STR$(ABS(Tz))+"hrs "+NS$+" Greenwich) On "+mydate$
END FUNCTION

Jim

VK7JH
MMedit

MMAndy
Regular Member

Joined: 16/07/2015
Location: United States
Posts: 91

Posted: 10:00pm 12 Sep 2015

Copy link to clipboard

Print this post

Quote I think you had the date order wrong.

Leaving out the commas will cause the wrong results.

WhiteWizzard
Guru

Joined: 05/04/2013
Location: United Kingdom
Posts: 2944

Posted: 10:26pm 12 Sep 2015

Copy link to clipboard

Print this post

If it helps, then to reduce 'input error' then modify:

. . . .
INPUT "Year, Month, Day";Y,M,D
IF Y<100 THEN Y = Y + 2000 ' allow for two digit years

to:

. . . .
getDate:
Y=-1:M=-1:D=-1 'reset variables to indicate missing Y,M,D entries
INPUT "Year, Month, Day";Y,M,D 'set variables; any missing will still =-1
IF y<0 OR M<0 OR D<0 THEN GOTO getDate 'check Y and M and D have been entered; if not then re-enter
IF Y<100 THEN Y = Y + 2000 'allow for two digit years

isochronic
Guru

Joined: 21/01/2012
Location: Australia
Posts: 689

Posted: 11:29pm 12 Sep 2015

Copy link to clipboard

Print this post

This is cut-and-paste descriptions of cos from the X32 math libraries doc:

cos
Description: Calculates the trigonometric cosine function of a double precision
floating-point value.
Include: <math.h>
Prototype: double cos (double x);
Argument: x value for which to return the cosine
Return Value: Returns the cosine of x in radians in the ranges of -1 to 1 inclusive.
Remarks: A domain error will occur if x is a NaN or infinity.

cosf
Description: Calculates the trigonometric cosine function of a single precision
floating-point value.
Include: <math.h>
Prototype: float cosf (float x);
Argument: x value for which to return the cosine
Return Value: Returns the cosine of x in radians in the ranges of -1 to 1 inclusive.
Remarks: A domain error will occur if x is a NaN or infinity.

So I think a C program using "cos" will use the dp function by default (promoting sp input), unless cosf is deliberately used instead. It is a strength, not a problem !

(ed - keeping track of what version of compiler does what with multiple libraries and chips and versions, is becoming an art form in itself )

Agreed, it is often not important but for solar concentrator work, a small aiming error is critical.
Edited by chronic 2015-09-14

MMAndy
Regular Member

Joined: 16/07/2015
Location: United States
Posts: 91

Posted: 02:49am 13 Sep 2015

Copy link to clipboard

Print this post

I see an application for the mighty MicroMite(+).
The above sunrise and sunset algorithm needs the year, month and day which can be obtained automatically from the MM+ external I2C RTC. (Date$)
For the time zone, which will change during the dates of daylight saving time (+- 1 hour) you can use the DST algorithm from Geoff Graham's computer controlled clock. He uses the date and day of week to compute the daylight saving time. You need to create a flag variable to indicate whether your in or out of daylight saving time. Now what you have is an algorithm which can be used to control indoor or outdoor lighting using TassyJim's sunrise and sunset time(s) algorithm.

Edited by MMAndy 2015-09-14

kiiid

Guru

Joined: 11/05/2013
Location: United Kingdom
Posts: 671

Posted: 04:59am 13 Sep 2015

Copy link to clipboard

Print this post

MMAndy said

...........

I don't think lighting normally depends strictly on hours, but rather on ambient luminosity instead. Hence for your particular application a simple photodiode based circuit with hysteresis, would do the same job more accurately from human's point of view :)

Edited by kiiid 2015-09-14

http://rittle.org

--------------

MMAndy
Regular Member

Joined: 16/07/2015
Location: United States
Posts: 91

Posted: 08:39am 13 Sep 2015

Copy link to clipboard

Print this post

I would use both a photodiode and the sunrise sunset algorithm.

As for an output, a SSR like a 3v input Crydom brick will drive any AC mains circuit.

BTW ... The sunrise sunset algorithm is awesome - TassyJim. Edited by MMAndy 2015-09-14

TassyJim

Guru

Joined: 07/08/2011
Location: Australia
Posts: 6283

Posted: 02:00pm 13 Sep 2015

Copy link to clipboard

Print this post

If you are looking a light levels, perhaps twilight is more appropriate.
There are 3 main definitions of twilight. Civil, Nautical and astronomical.

For Civil twilight at mu location we have:

Lat, Long (deg)? -41.125833,145.907778
Time zone (hrs)? -10
Year, Month, Day? 15,9,14

41deg 8min South, 145deg 54min West
(10hrs Ahead of Greenwich) On 14-Sep-2015

Sunrise at 06:21, azimuth 85.852
Sunset at 18:04, azimuth 273.903

Twilight:

Sunrise at 05:53, azimuth 90.3354
Sunset at 18:32, azimuth 269.386

All it takes is a small change to the function.
The extra variable passed to it is 0=sunset/rise, 1= civil, 2=nautical, 3=astronomical.

' Sunrise-Sunset
' TassyJim 13 Sep 2015
' Suitable for MMBasic V4.5 and above on Maximites and Micromites
'
' Latitude South and Longitude West are negative
' Time zones east of Greenwich (Australian) are negative
'
' Original program notes:
' This program by Roger W. Sinnott calculates the times of sunrise
' and sunset on any date, accurate to the minute within several
' centuries of the present. It correctly describes what happens in the
' arctic and antarctic regions, where the Sun may not rise or set on
' a given date. Enter north latitudes positive, west longitudes
' negative. For the time zone, enter the number of hours west of
' Greenwich (e.g., 5 for EST, 4 for EDT). The calculation is
' discussed in Sky & Telescope for August 1994, page 84.

IF MM.VER >4.05 THEN ' OPTION EXPLICIT introduced in V4.6
OPTION EXPLICIT
DIM Lat,Lon,Tz,D,M,Y
ENDIF

DO
INPUT "Lat, Long (deg)";Lat,Lon
INPUT "Time zone (hrs)";Tz
INPUT "Year, Month, Day";Y,M,D
IF Y<100 THEN Y = Y + 2000 ' allow for two digit years
PRINT ""
PRINT checkinput$(Lat,Lon,Tz,Y,M,D)
PRINT ""
PRINT suntimes$(Lat,Lon,Tz,Y,M,D,0)
PRINT "Twilight:"
PRINT suntimes$(Lat,Lon,Tz,Y,M,D,1)' 1 = Civil, 2=Nautical, 3 = Astronomical
PRINT ""

LOOP
END

FUNCTION suntimes$(Lat,Lon,Tz,Y,M,D,Tw)
' returns a string with sunrise and sunset times
LOCAL Ax(2), Dx(2), crlf$, Lonx, Dy
LOCAL K1,Z0,G,D1,F,J,S,A,J3,T,TT,T0,n,L,V,U,W,P,B
LOCAL A5,A0,D0,DA,DD,A2,D2,H0,H1,H2,L0,L2,C0,C,H7
LOCAL Z,M8,W8,V0,V2,E,T3,H3,M3,N7,D7,AZ,D5,Z1,V1, Za
LOCAL Sset$,Srise$,Msg1$,Msg2$,Msg3$,Msg4$
Dy=D
crlf$=CHR$(13)+CHR$(10)
K1=15*(PI/180)*1.0027379
Sset$="Sunset at "
Srise$="Sunrise at "
Msg1$="No sunrise this date"
Msg2$="No sunset this date"
Msg3$="Sun down all day"
Msg4$="Sun up all day"
IF Tw = 0 THEN ' sunrise/sunset
Za = 0.833
ELSE
Za = Tw*6 ' 1 = Civil, 2=Nautical, 3 = Astronomical
ENDIF

Lonx=Lon/360: Z0=Tz/24
G=1: IF Y<1583 THEN G=0
D1=INT(Dy): F=Dy-D1-.5
J=0-INT(7*(INT((M+9)/12)+Y)/4)
IF G<>0 THEN
S=SGN(M-9): A=ABS(M-9)
J3=INT(Y+S*INT(A/7))
J3=0-INT((INT(J3/100)+1)*3/4)
ENDIF
J=J+INT(275*M/9)+D1+G*J3
J=J+1721027+2*G+367*Y
IF F<0 THEN
F=F+1: J=J-1
ENDIF
T=(J-2451545)+F
TT=T/36525+1: ' TT = centuries from 1900
' LST at 0h zone time
T0=T/36525
S=24110.5+8640184.813*T0
S=S+86636.6*Z0+86400*Lonx
S=S/86400: S=S-INT(S)
T0=S*360*PI/180
T=T+Z0
' Get Sun's Position
FOR n = 1 TO 2
' Fundamental arguments
' (Van Flandern &
' Pulkkinen, 1979)
L=.779072+.00273790931*T
G=.993126+.0027377785*T
L=L-INT(L): G=G-INT(G)
L=L*PI*2: G=G*PI*2
V=.39785*SIN(L)
V=V-.01000*SIN(L-G)
V=V+.00333*SIN(L+G)
V=V-.00021*TT*SIN(L)
U=1-.03349*COS(G)
U=U-.00014*COS(2*L)
U=U+.00008*COS(L)
W=-.00010-.04129*SIN(2*L)
W=W+.03211*SIN(G)
W=W+.00104*SIN(2*L-G)
W=W-.00035*SIN(2*L+G)
W=W-.00008*TT*SIN(G)
' Compute Sun's RA and Dec
S=W/SQR(U-V*V)
A5=L+ATN(S/SQR(1-S*S))
S=V/SQR(U):D5=ATN(S/SQR(1-S*S))
'R5=1.00021*SQR(U) ' we don't use R5 anywhere
Ax(n)=A5: Dx(n)=D5
T=T+1
NEXT n
T=T-1
IF Ax(2)<Ax(1) THEN Ax(2)=Ax(2)+PI*2
Z1=(PI/180)*(90+Za): ' Zenith dist.
S=SIN(Lat*(PI/180)): C=COS(Lat*(PI/180))
Z=COS(Z1): M8=0: W8=0: PRINT
A0=Ax(1): D0=Dx(1)
DA=Ax(2)-Ax(1): DD=Dx(2)-Dx(1)
FOR C0=0 TO 23
P=(C0+1)/24
A2=Ax(1)+P*DA: D2=Dx(1)+P*DD
' Test an hour for an event
L0=T0+C0*K1: L2=L0+K1
H0=L0-A0: H2=L2-A2
H1=(H2+H0)/2: ' Hour angle,
D1=(D2+D0)/2: ' declination,
' at half hour
IF C0<=0 THEN
V0=S*SIN(D0)+C*COS(D0)*COS(H0)-Z
ENDIF
V2=S*SIN(D2)+C*COS(D2)*COS(H2)-Z
IF SGN(V0)<>SGN(V2) THEN
V1=S*SIN(D1)+C*COS(D1)*COS(H1)-Z
A=2*V2-4*V1+2*V0: B=4*V1-3*V0-V2
Dy=B*B-4*A*V0
IF Dy >= 0 THEN
Dy=SQR(Dy)
IF V0<0 AND V2>0 THEN suntimes$=suntimes$ + Srise$
IF V0<0 AND V2>0 THEN M8=1
IF V0>0 AND V2<0 THEN suntimes$=suntimes$ + Sset$
IF V0>0 AND V2<0 THEN W8=1
E=(0-B+Dy)/(2*A)
IF E>1 OR E<0 THEN E=(0-B-Dy)/(2*A)
T3=C0+E+1/120: ' Round off
H3=INT(T3): M3=INT((T3-H3)*60)
suntimes$=suntimes$ + ftime$(H3,M3)
H7=H0+E*(H2-H0)
N7=-COS(D1)*SIN(H7)
D7=C*SIN(D1)-S*COS(D1)*COS(H7)
AZ=ATN(N7/D7)/(PI/180)
IF D7<0 THEN AZ=AZ+180
IF AZ<0 THEN AZ=AZ+360
IF AZ>360 THEN AZ=AZ-360
suntimes$=suntimes$ + ", azimuth "+STR$(AZ)+crlf$
ENDIF
ENDIF
A0=A2: D0=D2: V0=V2
NEXT
IF M8=0 AND W8=0 THEN ' Special msg?
IF V2<0 THEN suntimes$=suntimes$ + Msg3$+crlf$
IF V2>0 THEN suntimes$=suntimes$ + Msg4$+crlf$
ELSE
IF M8=0 THEN suntimes$=suntimes$ + Msg1$+crlf$
IF W8=0 THEN suntimes$=suntimes$ + Msg2$+crlf$
ENDIF
END FUNCTION

FUNCTION ftime$(h,m) ' a time format method that works on most versions of MMBasic
ftime$=RIGHT$(STR$(h+100),2)+":"+RIGHT$(STR$(m+100),2)
END FUNCTION

FUNCTION checkinput$(Lat,Lon,Tz,Y,M,D)
LOCAL mydate$, mylat$, mylon$, NS$, de, mm, La, Lo
IF M>0 AND M<=12) THEN
mydate$=STR$(D)+"-"+MID$("xxJanFebMarAprMayJunJulAugSepOctNovDec",M*3,3)+"-"+STR$(Y)
ELSE
mydate$="Unknown date format"
ENDIF
IF D < 1 OR D >31 THEN
mydate$="Unknown date format"
ENDIF
IF Lat<0 THEN
NS$="South"
ELSE
NS$="North"
ENDIF
La = ABS(Lat)
de=INT(La): mm = INT((La - de)*60 + 0.5)
mylat$=STR$(de)+"deg "+STR$(mm)+"min "+NS$
IF Lon<0 THEN
NS$="East"
ELSE
NS$="West"
ENDIF
Lo = ABS(Lon)
de=INT(Lo): mm = INT((Lo - de)*60 + 0.5)
mylon$=STR$(de)+"deg "+STR$(mm)+"min "+NS$
IF Tz < 0 THEN
NS$="Ahead of"
ELSE
NS$="Behind"
ENDIF
checkinput$=mylat$+", "+mylon$+CHR$(13)+CHR$(10)
checkinput$=checkinput$+"("+STR$(ABS(Tz))+"hrs "+NS$+" Greenwich) On "+mydate$
END FUNCTION

Jim

VK7JH
MMedit

Print this page

To reply to this topic, you need to log in.