Search

Untitled

Here is an implementation of the Sunrise/Sunset Algorithm in VBA.

Official times of sunrise og sunset are calculated with variables depending on location given in latitude and longitude and corrected for difference between UCT/GMT and local time. Tested in Excel 2007.

Source of pseudocode algorithm: Almanac for Computers, 1990; Nautical Almanac Office; United States Naval Observatory; Washington, DC 20392

Inputs: day (integer), month (integer), year (double), latitude (double), longitude (double), localOffset (difference between local time and UCT/GMT), status (return time of official sun rise=0; return time of official sun set=1)

day, month and year are return values of standard Excel functions DAY, MONTH and YEAR applied to a cell formatted as an date.

My own location: 57.07 north and 9.94 east (eastern coordinates are positive, western negative) This year the sun will rise at 5:49:04 AM and set at 9:00:03 PM at my home address on my birthday.

Option Explicit

Public Const cPI As Double = 3.14159265358979

Public Function IntervalCorrection(a, b)     Dim Temp As Double     Temp = a     If Temp < 0 Then     Temp = Temp + b     End If     If Temp > b Then     Temp = Temp – b     End If     IntervalCorrection = Temp End Function

‘Implementation of trigonometric functions taking degrees as argument instead of radians Public Function cosDeg(degree As Double) cosDeg = Cos(degree * cPI / 180#) End Function Public Function sinDeg(degree As Double) sinDeg = Sin(degree * cPI / 180#) End Function Public Function tanDeg(degree As Double) tanDeg = Tan(degree * cPI / 180#) End Function ‘Implementation of trigonometric functions returning degrees instead of radians Public Function acosDeg(arg As Double) acosDeg = 180# * Application.WorksheetFunction.Acos(arg) / cPI End Function Public Function asinDeg(arg As Double) asinDeg = 180# * Application.WorksheetFunction.Asin(arg) / cPI End Function Public Function atanDeg(arg As Double) atanDeg = 180# * Atn(arg) / cPI End Function

Public Function suntimes(day As Integer, month As Integer, year As Integer, latitude As Double, longitude As Double, localOffset As Integer, state As Double) ‘The day of the year and associated intermediate variables Dim N As Double Dim N1 As Double Dim N2 As Double Dim N2A As Double Dim N3 As Double ‘Hour value and approximate time Dim lngHour As Double Dim t As Double ‘Suns mean anomaly Dim M As Double ‘Suns true longitude and intermediate variables Dim L As Double Dim L1 As Double Dim L2 As Double ‘Suns right ascension Dim RA As Double ‘Intermediate variables in calculation of ascension value Dim Lquadrant As Double Dim RAquadrant As Double ‘Suns declination Dim sinDec As Double Dim cosDec As Double ‘Suns local hour value and intermediate values Dim H As Double Dim cosH As Double Dim cosH1 As Double Dim cosH2 As Double Dim cosH3 As Double ‘Local mean time of rising Dim LT As Double ‘Coordinated Universal Time of rising Dim UT As Double ‘Definition of zenith – official 90.0 degrees, civil 96 50’ degrees, nautical 102 degrees, astronomical 108 degrees Dim zenith As Double zenith = 90#

‘Calculate the day of the year N1 = WorksheetFunction.Floor(275# * month / 9#, 1) N2 = WorksheetFunction.Floor((month + 9#) / 12#, 1) N2A = WorksheetFunction.Floor(year / 4#, 1) N3 = 1 + WorksheetFunction.Floor((year – 4# * N2A + 2#) / 3#, 1) N = N1 – (N2 * N3) + day – 30# Debug.Print “The value of variable N is: ” & N ‘Convert the longitude to hour value and calculate an approximate time lngHour = longitude / 15# t = N + ((6# – lngHour) / 24#) If (state > 0) Then t = N + ((18# – lngHour) / 24#) End If

‘Calculate the Sun’s mean anomaly M = (0.9856 * t) – 3.289 Debug.Print “The value of variable M is: ” & M ‘Calculate the Suns true longitude L1 = 1.916 * sinDeg(M) L2 = 0.02 * sinDeg(2# * M) L = M + L1 + L2 + 282.634 L = IntervalCorrection(L, 360)

‘Calculate the Suns right ascension RA = atanDeg(0.91764 * tanDeg(L)) RA = IntervalCorrection(RA, 360)

‘Right ascension needs to be in the same quadrant as L Lquadrant = (Application.WorksheetFunction.Floor(L / 90#, 1)) * 90# RAquadrant = (Application.WorksheetFunction.Floor(RA / 90#, 1)) * 90# RA = RA + (Lquadrant – RAquadrant)

‘Right ascension value needs to be converted into hours RA = RA / 15#

‘Calculate the Suns declination sinDec = 0.39782 * sinDeg(L) cosDec = cosDeg(asinDeg(sinDec))

‘Calculate the Suns local hour angle cosH1 = cosDeg(zenith) cosH2 = sinDec * sinDeg(latitude) cosH3 = cosDec * cosDeg(latitude) cosH = (cosH1 – cosH2) / cosH3

If (cosH > 1) Then     MsgBox (“The sun do not rise”)     Exit Function  End If

 If (cosH < -1) Then     MsgBox (“The sun do not set”)     Exit Function  End If

‘Finish calculating H and convert into hours H = (360 – acosDeg(cosH)) / 15# If (state > 0) Then H = acosDeg(cosH) / 15# End If

‘Calculate local mean time of rising LT = H + RA – (0.06571 * t) – 6.622 ‘Adjust back to UTC UT = LT – lngHour UT = IntervalCorrection(UT, 24) Debug.Print “The value of variable cPI is: ” & cPI Debug.Print “The value of variable UT is: ” & UT Debug.Print “The value of variable localOffset is: ” & localOffset ‘Return value suntimes = UT + localOffset End Function

0 views

Recent Posts

See All

©2020 by Danish Institute for Data Science. Proudly created with Wix.com