here is my last library for Solar calculations.
I use it for turn on/off my home garden lamps.
Enjoy,
Max
Code: Select all
'+-----------------------------------------------------------------------------+
'+ SunCalc_Library
'+
'+ Bytex (Italy)
'+ September 2008
'+
'+ Hardware configuration:
'+ Any PIC 18 family with 32KB ROM
'+ OSC: 8 ~ 20MHz
'+
'+ Features:
'+
'+ SunRise time calculation
'+ Noon time calculation
'+ SunSet time calculation
'+ Day of year calculation
'+ Day of week calculation
'+ Is Leap Year calculation
'+ Daylight Saving Time period, starting date, ending date calculation
'+
'+ web references from:
'+ http://www.srrb.noaa.gov/highlights/sunrise/sunrise.html
'+
'+-----------------------------------------------------------------------------+
module SunCalc_Library
'===============================================================================
' SC_ (SunCalc) GLOBAL INPUT VARIABLES
'===============================================================================
'+-------------------------------------+
'+ input parameters for local latitude +
'+ SC_in_LatDeg --> North +; South -; +
'+-------------------------------------+
Dim SC_in_LatDeg as Integer
Dim SC_in_LatMin As Integer
Dim SC_in_LatSec As Integer
Dim SC_in_LatSec100 As Integer
Dim SC_in_Latitude as float
'+--------------------------------------+
'+ input parameters for local longitude +
'+ SC_in_LonDeg --> West +; East -; +
'+--------------------------------------+
Dim SC_in_LonDeg As Integer
Dim SC_in_LonMin As Integer
Dim SC_in_LonSec As Integer
Dim SC_in_LonSec100 As Integer
Dim SC_in_Longitude as float
'+------------------------------------------------+
'+ input parameter for Greenwich Mean Time offset +
'+ SC_in_GMToffset --> West -; East +; +
'+------------------------------------------------+
Dim SC_in_GMToffset as integer
'+---------------------------------+
'+ input parameters for today date +
'+---------------------------------+
Dim SC_in_ThisDay As Integer
Dim SC_in_ThisMonth As Integer
Dim SC_in_ThisYear As Integer
'+--------------------------------------------------+
'+ input parameters for Daylight Saving Time period +
'+ For Europe starting period --> march +
'+ For Europe ending period ----> october +
'+--------------------------------------------------+
Dim SC_in_dst_SMonth as integer
Dim SC_in_dst_EMonth as integer
'===============================================================================
' SC_ (SunCalc) GLOBAL OUTPUT VARIABLES
'===============================================================================
'+---------------------------------+
'+ Output values for SunRise event +
'+---------------------------------+
Dim SC_out_RiseTime_Hours as byte
Dim SC_out_RiseTime_Minutes as byte
'+--------------------------------+
'+ Output values for SunSet event +
'+--------------------------------+
Dim SC_out_SetTime_Hours as byte
Dim SC_out_SetTime_Minutes as byte
'+-----------------------------------+
'+ Output values for Noon time event +
'+-----------------------------------+
Dim SC_out_NoonTime_Hours as byte
Dim SC_out_NoonTime_Minutes as byte
'+---------------------------------------------+
'+ String output for SunRise, NoonTime, SunSet +
'+ showing hh:mm +
'+---------------------------------------------+
Dim SC_out_StrBuffer as string[5]
'+----------------------------------------+
'+ Output value for numerical day of year +
'+ indicating 1..365 +
'+----------------------------------------+
Dim SC_out_DayOfYear as integer
'+-------------------------------------+
'+ String output value for day of week +
'+ 1=monday...7=sunday +
'+-------------------------------------+
Dim SC_out_StrDayOfWeek as string[3]
'+------------------------------+
'+ Output value for day of week +
'+ 1=monday...7=sunday +
'+------------------------------+
Dim SC_out_IntDayOfWeek as integer
'+-----------------------------------+
'+ Output value for leap year result +
'+ yes=true; no=false +
'+-----------------------------------+
Dim SC_out_IsLeapYear as byte
'+-----------------------------------------------------+
'+ Output value for Daylight Saving Time period result +
'+ yes=true; no=false +
'+-----------------------------------------------------+
Dim SC_out_DaylightSavingTime as byte
'+----------------------------------------+
'+ Output values for Daylight Saving Time +
'+ starting date +
'+----------------------------------------+
Dim SC_out_DST_StartDay as integer
Dim SC_out_DST_StartMonth as integer
Dim SC_out_DST_StartYear as integer
'+----------------------------------------+
'+ Output values for Daylight Saving Time +
'+ ending date +
'+----------------------------------------+
Dim SC_out_DST_EndDay as integer
Dim SC_out_DST_EndMonth as integer
Dim SC_out_DST_EndYear as integer
implements
'//****************************************************************************/
'//* Name: GetLocalLatitude */
'//* Type: Function */
'//* Purpose: */
'//* Arguments: */
'//* */
'//****************************************************************************/
Sub Function GetLocalLatitude(Dim SC_in_LatDeg As Integer,
Dim SC_in_LatMin As Integer,
Dim SC_in_LatSec As Integer,
Dim SC_in_LatSec100 As Integer) As float
Dim w As float
Dim p As float
p = float(SC_in_LatSec) + float(SC_in_LatSec100) / 100.0
If SC_in_LatDeg > 0 Then
w = float(SC_in_LatDeg) / 1.0
w = w + float(SC_in_LatMin) / 60.0
w = w + p / 3600.0
Else
w = float(SC_in_LatDeg) / 1.0
w = w - float(SC_in_LatMin) / 60.0
w = w - p / 3600.0
End If
If w < -89.0 Then
w = -89.0
End If
If w > 89.0 Then
w = 89.0
End If
Result = w
End sub
'//****************************************************************************/
'//* Name: GetSC_in_Longitude */
'//* Type: Function */
'//* Purpose: */
'//* Arguments: */
'//* */
'//****************************************************************************/
Sub Function GetLocalLongitude(Dim SC_in_LonDeg As Integer,
Dim SC_in_LonMin As Integer,
Dim SC_in_LonSec As Integer,
Dim SC_in_LonSec100 As Integer) as Float
Dim w as Float
Dim p as Float
p = float(SC_in_LonSec) + float(SC_in_LonSec100) / 100.0
If SC_in_LonDeg > 0 Then
w = float(SC_in_LonDeg) / 1.0
w = w + float(SC_in_LonMin) / 60.0
w = w + p / 3600.0
Else
w = float(SC_in_LonDeg) / 1.0
w = w - float(SC_in_LonMin) / 60.0
w = w - p / 3600.0
End If
Result = w
End Sub
'//****************************************************************************/
'//* Name: isLeapYear */
'//* Type: Function */
'//* Purpose: */
'//* Arguments: */
'//* */
'//****************************************************************************/
Sub Function isLeapYear(Dim yr As Integer) As byte
If ((yr Mod 4 = 0) And (yr Mod 100 <> 0)) Or (yr Mod 400 = 0) Then
Result = True
Else
Result = False
End If
End Sub
'=====================================================================
' Convert radian angle to degrees
'=====================================================================
Sub Function RadToDeg(Dim RadAngle as Float) as Float
Result = (180.0 * RadAngle / 3.14159265358979)
End Sub
'=====================================================================
' Convert degree angle to radians
'=====================================================================
Sub Function DegToRad(Dim DegAngle as Float) as Float
Result = (DegAngle * 3.14159265358979 / 180.0)
End Sub
'//****************************************************************************/
'//* Name: calcDayOfYear */
'//* Type: Function */
'//* Purpose: Finds numerical day-of-year from mn, day and lp year info */
'//* Arguments: */
'//* month: January = 1 */
'//* day : 1 - 31 */
'//* lpyr : 1 if leap year, 0 if not */
'//* Return value: */
'//* The numerical day of year */
'//****************************************************************************/
Sub Function calcDayOfYear(Dim mMonth As Integer,
Dim mday As Integer,
Dim lpyr As byte) As integer
Dim k As Integer
If lpyr = True Then
k = 1
Else
k = 2
End If
Result = integer(floor((275.0 * float(mMonth)) / 9.0) -
float(k) * floor((float(mMonth) + 9.0) /
12.0) + float(mday) - 30.0)
End Sub
'//***********************************************************************/
'//* Name: calcDayOfWeek */
'//* Type: Function */
'//* Purpose: Derives weekday from Julian Day */
'//* Arguments: */
'//* juld : Julian Day */
'//* Return value: */
'//* String containing name of weekday */
'//***********************************************************************/
Sub Procedure calcDayOfWeek(Dim JulD as Float,
Dim ByRef StrDayOfWeek as string[3],
Dim Byref IntDayOfWeek as integer)
Dim a1,a2 as float
Dim ba as float
Dim bc as integer
StrDayOfWeek = " "
IntDayOfWeek = 0
a1 = JulD + 1.5
a2 = a1 / 7.0
ba = a1 - floor(a2) * 7.0
bc = integer(ba)
Select Case bc
Case 0
StrDayOfWeek = "Sun"
Case 1
StrDayOfWeek = "Mon"
Case 2
StrDayOfWeek = "Tue"
Case 3
StrDayOfWeek = "Wed"
Case 4
StrDayOfWeek = "Thu"
Case 5
StrDayOfWeek = "Fri"
Case 6
StrDayOfWeek = "Sat"
End Select
if bc > 0 then
IntDayOfWeek = bc
else
IntDayOfWeek = 7
end if
End Sub
'//*** [1] ***************************************************************/
'//* Name: calcJD */
'//* Type: Function */
'//* Purpose: Julian day from calendar day */
'//* Arguments: */
'//* year : 4 digit year */
'//* month: January = 1 */
'//* day : 1 - 31 */
'//* Return value: */
'//* The Julian day corresponding to the date */
'//* Note: */
'//* Number is returned for start of day. Fractional days should be */
'//* added later. */
'//***********************************************************************/
Sub Function CalcJD(dim aYear As Integer,
dim aMonth As Integer,
dim aDay As Integer) As float
'+ Web references ---------------------------+
'+ +
'+ http://en.wikipedia.org/wiki/Julian_day +
'+ +
'+-------------------------------------------+
Dim cMonth As integer
Dim cYear As integer
Dim cDay As integer
Dim A As integer
Dim B As integer
If aMonth <= 2 Then
cYear = aYear - 1
cMonth = aMonth + 12
Else
cYear = aYear
cMonth = aMonth
End If
cDay = aDay
A = floor(cYear / 100.0)
B = 2.0 - A + floor(A / 4.0)
Result = floor(365.25 * (cYear + 4716)) +
floor(30.6001 * (cMonth + 1)) +
cDay + B - 1524.5
End Sub
'//*** [2] **************************************************************/
'//* Name: calcTimeJulianCent */
'//* Type: Function */
'//* Purpose: convert Julian Day to centuries since J2000.0. */
'//* Arguments: */
'//* jd : the Julian Day to convert */
'//* Return value: */
'//* the value corresponding to the Julian Day */
'//**********************************************************************/
Sub Function CalcTimeJulianCent(Dim JD as Float) as Float
Result = (JD - 2451545.0) / 36525.0
End Sub
'//**********************************************************************/
'//* Name: calcJDFromJulianCent */
'//* Type: Function */
'//* Purpose: convert centuries since J2000.0 to Julian Day. */
'//* Arguments: */
'//* t : number of Julian centuries since J2000.0 */
'//* Return value: */
'//* the Julian Day corresponding to the t value */
'//**********************************************************************/
Sub Function CalcJDFromJulianCent(Dim t as Float) as Float
Result = t * 36525.0 + 2451545.0
End Sub
'//**********************************************************************/
'//* Name: calGeomMeanLongSun */
'//* Type: Function */
'//* Purpose: calculate the Geometric Mean SC_in_Longitude of the Sun */
'//* Arguments: */
'//* t : number of Julian centuries since J2000.0 */
'//* Return value: */
'//* the Geometric Mean SC_in_Longitude of the Sun in degrees */
'//**********************************************************************/
Sub Function CalcGeomMeanLongSun(Dim t as Float) as Float
Dim L0 as Float
L0 = 280.46646 + t * (36000.76983 + 0.0003032 * t)
While (L0 > 360.0)
L0 = L0 - 360.0
Wend
While (L0 < 0.0)
L0 = L0 + 360.0
Wend
Result = L0
End Sub
'//**********************************************************************/
'//* Name: calGeomAnomalySun */
'//* Type: Function */
'//* Purpose: calculate the Geometric Mean Anomaly of the Sun */
'//* Arguments: */
'//* t : number of Julian centuries since J2000.0 */
'//* Return value: */
'//* the Geometric Mean Anomaly of the Sun in degrees */
'//**********************************************************************/
Sub Function CalcGeomMeanAnomalySun(Dim t as Float) as Float
Result = 357.52911 + t * (35999.05029 - 0.0001537 * t)
End Sub
'//**********************************************************************/
'//* Name: calcEccentricityEarthOrbit */
'//* Type: Function */
'//* Purpose: calculate the eccentricity of earth's orbit */
'//* Arguments: */
'//* t : number of Julian centuries since J2000.0 */
'//* Return value: */
'//* the unitless eccentricity */
'//**********************************************************************/
Sub Function CalcEccentricityEarthOrbit(Dim t as Float) as Float
Result = 0.016708634 - t * (0.000042037 + 0.0000001267 * t)
End Sub
'//**********************************************************************/
'//* Name: calcMeanObliquityOfEcliptic */
'//* Type: Function */
'//* Purpose: calculate the mean obliquity of the ecliptic */
'//* Arguments: */
'//* t : number of Julian centuries since J2000.0 */
'//* Return value: */
'//* mean obliquity in degrees */
'//**********************************************************************/
Sub Function CalcMeanObliquityOfEcliptic(Dim t as Float) as Float
Dim seconds as Float
seconds = 21.448 - t * (46.815 + t * (0.00059 - t * (0.001813)))
Result = 23.0 + (26.0 + (seconds / 60.0)) / 60.0
End Sub
'//**********************************************************************/
'//* Name: calcObliquityCorrection */
'//* Type: Function */
'//* Purpose: calculate the corrected obliquity of the ecliptic */
'//* Arguments: */
'//* t : number of Julian centuries since J2000.0 */
'//* Return value: */
'//* corrected obliquity in degrees */
'//**********************************************************************/
Sub Function CalcObliquityCorrection(Dim t as Float) as Float
Dim omega as Float
omega = 125.04 - 1934.136 * t
Result = CalcMeanObliquityOfEcliptic(t) + 0.00256 * Cos(DegToRad(omega))
End Sub
'//**********************************************************************/
'//* Name: calcEquationOfTime */
'//* Type: Function */
'//* Purpose: calculate the difference between true solar time and mean */
'//* solar time */
'//* Arguments: */
'//* t : number of Julian centuries since J2000.0 */
'//* Return value: */
'//* equation of time in minutes of time */
'//**********************************************************************/
Sub Function CalcEquationOfTime(Dim t as Float) as Float
Dim Epsilon as Float
Dim L0 as Float
Dim e as Float
Dim m as Float
Dim y as Float
Dim Sin2l0 as Float
Dim Sinm as Float
Dim Cos2l0 as Float
Dim Sin4l0 as Float
Dim Etime as Float
Epsilon = CalcObliquityCorrection(t)
L0 = CalcGeomMeanLongSun(t)
e = CalcEccentricityEarthOrbit(t)
m = CalcGeomMeanAnomalySun(t)
y = Tan(DegToRad(Epsilon) / 2.0)
y = y * y
Sin2l0 = Sin(2.0 * DegToRad(L0))
Sinm = Sin(DegToRad(m))
Cos2l0 = Cos(2.0 * DegToRad(L0))
Sin4l0 = Sin(4.0 * DegToRad(L0))
Etime = y * Sin2l0 - 2.0 * e * Sinm +
4.0 * e * y * Sinm *
Cos2l0 - 0.5 * y * y *
Sin4l0 - 1.25 * e * e *
Sin(2.0 * DegToRad(m))
Result = RadToDeg(Etime) * 4.0
End Sub
'//**********************************************************************/
'//* Name: calcSolNoonUTC */
'//* Type: Function */
'//* Purpose: calculate the Universal Coordinated Time (UTC) of solar */
'//* noon for the given day at the given location on earth */
'//* Arguments: */
'//* t : number of Julian centuries since J2000.0 */
'//* SC_in_Longitude : SC_in_Longitude of observer in degrees */
'//* Return value: */
'//* time in minutes from zero Z */
'//**********************************************************************/
Sub Function CalcSolNoonUTC(Dim t as Float,
Dim SC_in_Longitude as Float) as Float
Dim tnoon as Float
Dim eqTime as Float
Dim solNoonUTC as Float
Dim NewT as Float
tnoon = CalcTimeJulianCent(CalcJDFromJulianCent(t) + SC_in_Longitude / 360.0)
eqTime = CalcEquationOfTime(tnoon)
solNoonUTC = 720.0 + (SC_in_Longitude * 4.0) - eqTime '// min
NewT = CalcTimeJulianCent(CalcJDFromJulianCent(t) -
0.5 + solNoonUTC / 1440.0)
eqTime = CalcEquationOfTime(NewT)
solNoonUTC = 720.0 + (SC_in_Longitude * 4.0) - eqTime ' // min
Result = solNoonUTC
End Sub
'//**********************************************************************/
'//* Name: calcSunEqOfCenter */
'//* Type: Function */
'//* Purpose: calculate the equation of center for the sun */
'//* Arguments: */
'//* t : number of Julian centuries since J2000.0 */
'//* Return value: */
'//* in degrees */
'//**********************************************************************/
Sub Function CalcSunEqOfCenter(Dim t as Float) as Float
Dim mRad as Float
Dim Sin1m as Float
Dim Sin2m as Float
Dim Sin3m as Float
Dim C as Float
mRad = DegToRad(CalcGeomMeanAnomalySun(t))
Sin1m = Sin(mRad)
Sin2m = Sin(mRad + mRad)
Sin3m = Sin(mRad + mRad + mRad)
C = Sin1m * (1.914602 - t * (0.004817 + 0.000014 * t)) +
Sin2m * (0.019993 - 0.000101 * t) + Sin3m * 0.000289
Result = C ' in degrees
End Sub
'//**********************************************************************/
'//* Name: calcSunTrueLong */
'//* Type: Function */
'//* Purpose: calculate the true SC_in_Longitude of the sun */
'//* Arguments: */
'//* t : number of Julian centuries since J2000.0 */
'//* Return value: */
'//* sun's true SC_in_Longitude in degrees */
'//**********************************************************************/
Sub Function CalcSunTrueLong(Dim t as Float) as Float
Result = CalcGeomMeanLongSun(t) + CalcSunEqOfCenter(t)
End Sub
'//**********************************************************************/
'//* Name: calcSunTrueAnomaly */
'//* Type: Function */
'//* Purpose: calculate the true anamoly of the sun */
'//* Arguments: */
'//* t : number of Julian centuries since J2000.0 */
'//* Return value: */
'//* sun's true anamoly in degrees */
'//**********************************************************************/
Sub Function CalcSunTrueAnomaly(Dim t as Float) as Float
Result = CalcGeomMeanAnomalySun(t) + CalcSunEqOfCenter(t)
End Sub
'//**********************************************************************/
'//* Name: calcSunRadVector */
'//* Type: Function */
'//* Purpose: calculate the distance to the sun in AU */
'//* Arguments: */
'//* t : number of Julian centuries since J2000.0 */
'//* Return value: */
'//* sun radius vector in AUs */
'//**********************************************************************/
Sub Function CalcSunRadVector(Dim t as Float) as Float
Dim v as Float
Dim ee as Float
v = CalcSunTrueAnomaly(t)
ee = CalcEccentricityEarthOrbit(t)
Result = (1.000001018 * (1.0 - ee * ee)) / (1.0 + ee * Cos(DegToRad(v)))
End Sub
'//**********************************************************************/
'//* Name: calcSunApparentLong */
'//* Type: Function */
'//* Purpose: calculate the apparent SC_in_Longitude of the sun */
'//* Arguments: */
'//* t : number of Julian centuries since J2000.0 */
'//* Return value: */
'//* sun's apparent SC_in_Longitude in degrees */
'//**********************************************************************/
Sub Function CalcSunApparentLong(Dim t as Float) as Float
Result = CalcSunTrueLong(t) - 0.00569 -
0.00478 * Sin(DegToRad(125.04 - 1934.136 * t))
End Sub
'//**********************************************************************/
'//* Name: calcSunRtAscension */
'//* Type: Function */
'//* Purpose: calculate the right ascension of the sun */
'//* Arguments: */
'//* t : number of Julian centuries since J2000.0 */
'//* Return value: */
'//* sun's right ascension in degrees */
'//**********************************************************************/
Sub Function CalcSunRtAscension(Dim t as Float) as Float
Dim e as Float
Dim lambda as Float
Dim tananum as Float
Dim tanadenom as Float
e = CalcObliquityCorrection(t)
lambda = CalcSunApparentLong(t)
tananum = (Cos(DegToRad(e)) * Sin(DegToRad(lambda)))
tanadenom = (Cos(DegToRad(lambda)))
Result = RadToDeg(Atan2(tananum, tanadenom))
End Sub
'//**********************************************************************/
'//* Name: calcSunDeclination */
'//* Type: Function */
'//* Purpose: calculate the declination of the sun */
'//* Arguments: */
'//* t : number of Julian centuries since J2000.0 */
'//* Return value: */
'//* sun's declination in degrees */
'//**********************************************************************/
Sub Function CalcSunDeclination(Dim t as Float) as Float
Dim e as Float
Dim lambda as Float
Dim sint as Float
e = CalcObliquityCorrection(t)
lambda = CalcSunApparentLong(t)
sint = Sin(DegToRad(e)) * Sin(DegToRad(lambda))
Result = RadToDeg(Asin(sint))
End Sub
'//**********************************************************************/
'//* Name: calcHourAngleSunrise */
'//* Type: Function */
'//* Purpose: calculate the hour angle of the sun at sunrise for the */
'//* SC_in_Latitude */
'//* Arguments: */
'//* lat : SC_in_Latitude of observer in degrees */
'//* solarDec : declination angle of sun in degrees */
'//* Return value: */
'//* hour angle of sunrise in radians */
'//**********************************************************************/
Sub Function CalcHourAngleSunrise(Dim lat as Float,
Dim SolarDec as Float) as Float
Dim latrad as Float
Dim sdrad as Float
Dim HAarg as Float
latrad = DegToRad(lat)
sdrad = DegToRad(SolarDec)
HAarg = (Cos(DegToRad(90.833)) / (Cos(latrad) * Cos(sdrad)) -
Tan(latrad) * Tan(sdrad))
Result = (Acos(Cos(DegToRad(90.833)) / (Cos(latrad) * Cos(sdrad)) -
Tan(latrad) * Tan(sdrad)))
End Sub
'//**********************************************************************/
'//* Name: calcHourAngleSunset */
'//* Type: Function */
'//* Purpose: calculate the hour angle of the sun at sunset for the */
'//* SC_in_Latitude */
'//* Arguments: */
'//* lat : SC_in_Latitude of observer in degrees */
'//* solarDec : declination angle of sun in degrees */
'//* Return value: */
'//* hour angle of sunset in radians */
'//**********************************************************************/
Sub Function CalcHourAngleSunset(Dim lat as Float,
Dim SolarDec as Float) as Float
Dim latrad as Float
Dim sdrad as Float
Dim HAarg as Float
Dim HA as Float
latrad = DegToRad(lat)
sdrad = DegToRad(SolarDec)
HAarg = (Cos(DegToRad(90.833)) / (Cos(latrad) * Cos(sdrad)) -
Tan(latrad) * Tan(sdrad))
HA = (Acos(Cos(DegToRad(90.833)) / (Cos(latrad) * Cos(sdrad)) -
Tan(latrad) * Tan(sdrad)))
Result = -HA
End Sub
'//*** [3] ********************************************************************/
'//* Name: calcSunriseUTC */
'//* Type: Function */
'//* Purpose: calculate the Universal Coordinated Time (UTC) of sunrise */
'//* for the given day at the given location on earth */
'//* Arguments: */
'//* JD : julian day */
'//* SC_in_Latitude : SC_in_Latitude of observer in degrees */
'//* SC_in_Longitude : SC_in_Longitude of observer in degrees */
'//* Return value: */
'//* time in minutes from zero Z */
'//****************************************************************************/
Sub Function CalcSunRiseUTC(Dim JD as Float,
Dim SC_in_Latitude as Float,
Dim SC_in_Longitude as Float) as float
Dim t as Float
Dim noonmin as Float
Dim tnoon as Float
Dim eqTime as Float
Dim SolarDec as Float
Dim HourAngle as Float
Dim Delta as Float
Dim TimeDiff as Float
Dim TimeUTC as Float
Dim NewT as Float
t = CalcTimeJulianCent(JD)
noonmin = CalcSolNoonUTC(t, SC_in_Longitude)
tnoon = CalcTimeJulianCent(JD + noonmin / 1440.0)
eqTime = CalcEquationOfTime(tnoon)
SolarDec = CalcSunDeclination(tnoon)
HourAngle = CalcHourAngleSunrise(SC_in_Latitude, SolarDec)
Delta = SC_in_Longitude - RadToDeg(HourAngle)
TimeDiff = 4.0 * Delta '// in minutes of time
TimeUTC = 720.0 + TimeDiff - eqTime '// in minutes
NewT = CalcTimeJulianCent(CalcJDFromJulianCent(t) + TimeUTC / 1440.0)
eqTime = CalcEquationOfTime(NewT)
SolarDec = CalcSunDeclination(NewT)
HourAngle = CalcHourAngleSunrise(SC_in_Latitude, SolarDec)
Delta = SC_in_Longitude - RadToDeg(HourAngle)
TimeDiff = 4.0 * Delta
TimeUTC = 720.0 + TimeDiff - eqTime '// in minutes
Result = TimeUTC
End Sub
'//**********************************************************************/
'//* Name: calcSunsetUTC */
'//* Type: Function */
'//* Purpose: calculate the Universal Coordinated Time (UTC) of sunset */
'//* for the given day at the given location on earth */
'//* Arguments: */
'//* JD : julian day */
'//* SC_in_Latitude : SC_in_Latitude of observer in degrees */
'//* SC_in_Longitude : SC_in_Longitude of observer in degrees */
'//* Return value: */
'//* time in minutes from zero Z */
'//**********************************************************************/
Sub Function CalcSunSetUTC(Dim JD as Float,
Dim SC_in_Latitude as Float,
Dim SC_in_Longitude as Float) as Float
Dim t as Float
Dim noonmin as Float
Dim tnoon as Float
Dim eqTime as Float
Dim SolarDec as Float
Dim HourAngle as Float
Dim Delta as Float
Dim TimeDiff as Float
Dim TimeUTC as Float
Dim NewT as Float
t = CalcTimeJulianCent(JD)
noonmin = CalcSolNoonUTC(t, SC_in_Longitude)
tnoon = CalcTimeJulianCent(JD + noonmin / 1440.0)
eqTime = CalcEquationOfTime(tnoon)
SolarDec = CalcSunDeclination(tnoon)
HourAngle = CalcHourAngleSunset(SC_in_Latitude, SolarDec)
Delta = SC_in_Longitude - RadToDeg(HourAngle)
TimeDiff = 4.0 * Delta
TimeUTC = 720.0 + TimeDiff - eqTime
NewT = CalcTimeJulianCent(CalcJDFromJulianCent(t) + TimeUTC / 1440.0)
eqTime = CalcEquationOfTime(NewT)
SolarDec = CalcSunDeclination(NewT)
HourAngle = CalcHourAngleSunset(SC_in_Latitude, SolarDec)
Delta = SC_in_Longitude - RadToDeg(HourAngle)
TimeDiff = 4.0 * Delta
TimeUTC = 720.0 + TimeDiff - eqTime '// in minutes
Result = TimeUTC
End Sub
Sub Procedure TimeToString(Dim fTime as Float,
Dim GMT As Integer,
Dim Byref Buffer as string[5],
Dim ByRef ThoursOut as byte,
Dim Byref TminutesOut as byte)
Dim fHour As float
Dim fMin As float
Dim fSec As float
dim sHour as string[17]
dim sMin as string[17]
dim n as byte
dim k as byte
Buffer = " "
fSec = floor((fTime - floor(fTime)) * 60.0)
fSec = fSec + floor(fTime) * 60.0
fHour = floor(fSec / 3600.0)
fMin = floor((fSec / 3600.0 - fHour) * 60.0)
fSec = fSec - fHour * 3600.0 - fMin * 60.0
If fSec >= 30.0 Then
fMin = fMin + 1.0
End If
if SC_out_DaylightSavingTime then
fHour = fHour + float(GMT) + 1.0
else
fHour = fHour + float(GMT)
end if
floatToStr(fHour,sHour)
floatToStr(fMin,sMin)
k = 0
for n = 0 to 17
if sHour[n] = "." then
k = n
end if
next n
if k = 1 then
Buffer[0] = "0"
Buffer[1] = sHour[0]
end if
if k = 2 then
Buffer[0] = sHour[0]
Buffer[1] = sHour[1]
end if
Buffer[2] = ":"
k = 0
for n = 0 to 17
if sMin[n] = "." then
k = n
end if
next n
if k = 1 then
Buffer[3] = "0"
Buffer[4] = sMin[0]
end if
if k = 2 then
Buffer[3] = sMin[0]
Buffer[4] = sMin[1]
end if
ThoursOut = fHour
TminutesOut = fMin
End Sub
'//**********************************************************************/
'//* Name: GetLocalSunRiseTime */
'//* Type: Procedure */
'//* Purpose: calculate the local Sun Rise time */
'//* */
'//* Arguments: */
'//* JD : julian day */
'//* SC_in_Latitude : SC_in_Latitude of observer in degrees */
'//* SC_in_Longitude : SC_in_Longitude of observer in degrees */
'//* Return value: */
'//* time in minutes from zero Z */
'//**********************************************************************/
Sub Procedure GetLocalSunRiseTime(Dim SC_in_Latitude as Float,
Dim SC_in_Longitude as Float,
Dim mYear As Integer,
Dim mMonth As Integer,
Dim mDay As Integer,
Dim SC_in_GMToffset as integer,
Dim ByRef StrOut as string[5],
Dim ByRef ThoursOut as byte,
Dim Byref TminutesOut as byte)
Dim JD as Float
Dim t as Float
Dim RiseTimeGMT as Float
JD = CalcJD(mYear, mMonth, mDay)
t = CalcTimeJulianCent(JD)
' Apparent SunRise GMT time
RiseTimeGMT = CalcSunRiseUTC(JD, SC_in_Latitude, SC_in_Longitude)
TimeToString (RiseTimeGMT, SC_in_GMToffset , StrOut, ThoursOut, TminutesOut)
End Sub
Sub Procedure GetLocalSunSetTime(Dim SC_in_Latitude as Float,
Dim SC_in_Longitude as Float,
Dim mYear As Integer,
Dim mMonth As Integer,
Dim mDay As Integer,
Dim SC_in_GMToffset as integer,
Dim ByRef StrOut as string[5],
Dim ByRef ThoursOut as byte,
Dim Byref TminutesOut as byte)
Dim JD as Float
Dim t as Float
Dim SetTimeGMT as Float
JD = CalcJD(mYear, mMonth, mDay)
t = CalcTimeJulianCent(JD)
' Apparent SunSet GMT time
SetTimeGMT = CalcSunSetUTC(JD, SC_in_Latitude, SC_in_Longitude)
TimeToString (SetTimeGMT, SC_in_GMToffset , StrOut, ThoursOut, TminutesOut)
End Sub
Sub Procedure GetLocalNoonTime(Dim SC_in_Longitude as Float,
Dim mYear As Integer,
Dim mMonth As Integer,
Dim mDay As Integer,
Dim SC_in_GMToffset as integer,
Dim ByRef StrOut as string[5],
Dim ByRef ThoursOut as byte,
Dim Byref TminutesOut as byte)
Dim JD as Float
Dim t as Float
Dim SolNoonGMT as Float
JD = CalcJD(mYear, mMonth, mDay)
t = CalcTimeJulianCent(JD)
' Calculate solar noon for this date
SolNoonGMT = CalcSolNoonUTC(t, SC_in_Longitude)
TimeToString (SolNoonGMT, SC_in_GMToffset , StrOut, ThoursOut, TminutesOut)
End Sub
Sub Procedure GetMoreCalendarData(Dim mYear As Integer,
Dim mMonth As Integer,
Dim mDay As Integer,
Dim Byref mIsLeapYear as byte,
Dim Byref DayOfYear as integer,
Dim Byref StrDayOfWeek as string[3],
Dim Byref IntDayOfWeek as integer,
Dim Byref SC_out_DaylightSavingTime as byte)
Dim Jd as float
Dim i as integer
Dim LastSunday as integer
Dim sTemp as string[3]
Dim iTemp as integer
'+----------------------------------------------+
'+ Leap year calculation +
'+----------------------------------------------+
mIsLeapYear = isLeapYear(mYear)
'+----------------------------------------------+
'+ Numerical day of year calculation +
'+----------------------------------------------+
DayOfYear = calcDayOfYear(mMonth,mDay,mIsLeapYear)
'+----------------------------------------------+
'+ Numerical and string day of week calculation +
'+----------------------------------------------+
jd = CalcJD(mYear,mMonth,mDay)
calcDayOfWeek(jd,StrDayOfWeek,IntDayOfWeek)
'+----------------------------------------------+
'+ Daylight Saving Time (DST) calculation +
'+----------------------------------------------+
' remember to put the right ending day
' for the for-next loop.
'------------------------------------------
for i = 1 to 31 ' march ending day
iTemp = 0
jd = 0
jd = CalcJD(mYear,SC_in_dst_SMonth,i)
calcDayOfWeek(jd,sTemp,iTemp)
if iTemp = 7 then ' Looking for Sunday day
LastSunday = i
end if
next i
SC_out_DST_StartDay = LastSunday
SC_out_DST_StartMonth = SC_in_dst_SMonth
SC_out_DST_StartYear = mYear
' remember to put the right ending day
' for the for-next loop.
'------------------------------------------
for i = 1 to 31 ' october ending day
iTemp = 0
jd = 0
jd = CalcJD(mYear,SC_in_dst_EMonth,i)
calcDayOfWeek(jd,sTemp,iTemp)
if iTemp = 7 then ' Looking for Sunday day
LastSunday = i
end if
next i
SC_out_DST_EndDay = LastSunday
SC_out_DST_EndMonth = SC_in_dst_EMonth
SC_out_DST_EndYear = mYear
SC_out_DaylightSavingTime = false
if (CalcJD(mYear,mMonth,mDay) >=
CalcJD(SC_out_DST_StartYear,
SC_out_DST_StartMonth,
SC_out_DST_StartDay)) then
if (CalcJD(mYear,mMonth,mDay) <=
CalcJD(SC_out_DST_EndYear,
SC_out_DST_EndMonth,
SC_out_DST_EndDay)) then
'-----------------
SC_out_DaylightSavingTime = true
end if
end if
End Sub
'****************************** End Module *************************************
end.
Code: Select all
program SunRiseSunSet
include "SunCalc_Library"
Symbol LED7 = PORTB.7
Symbol LED6 = PORTB.6
Sub Procedure ConfigureLocalData
'==================================================
' Local latitude, longitude
' 45 54 19.48 , 13 18 36.26
'==================================================
SC_in_LatDeg = 45 ' North +; South -;
SC_in_LatMin = 54
SC_in_LatSec = 19
SC_in_LatSec100 = 48
SC_in_LonDeg = -13 ' West +; East -;
SC_in_LonMin = 18
SC_in_LonSec = 36
SC_in_LonSec100 = 26
SC_in_GMToffset = 1 ' West -; East +;
SC_in_Latitude = GetLocalLatitude(SC_in_LatDeg,
SC_in_LatMin,
SC_in_LatSec,
SC_in_LatSec100)
SC_in_Longitude = GetLocalLongitude(SC_in_LonDeg,
SC_in_LonMin,
SC_in_LonSec,
SC_in_LonSec100)
' Daylight Saving Time parameters
'-----------------------------------------------------------
SC_in_dst_SMonth = 3 ' Start period: march for Europe
SC_in_dst_EMonth = 10 ' End period: october for Europe
End Sub
main:
ConfigureLocalData
TRISB = 0
PORTB = 0
' Today date
'----------------------
SC_in_ThisYear = 2008
SC_in_ThisMonth = 9
SC_in_ThisDay = 29
while true
' Get more useful information from calendar
'--------------------------------------------------
GetMoreCalendarData(SC_in_ThisYear,
SC_in_ThisMonth,
SC_in_ThisDay,
SC_out_IsLeapYear,
SC_out_DayOfYear,
SC_out_StrDayOfWeek,
SC_out_IntDayOfWeek,
SC_out_DaylightSavingTime)
if SC_out_DaylightSavingTime then
LED7 = true
ELSE
LED7 = false
end if
GetLocalSunRiseTime (SC_in_Latitude,
SC_in_Longitude,
SC_in_ThisYear,
SC_in_ThisMonth,
SC_in_ThisDay,
SC_in_GMToffset,
SC_out_StrBuffer,
SC_out_RiseTime_Hours,
SC_out_RiseTime_Minutes)
GetLocalNoonTime (SC_in_Longitude,
SC_in_ThisYear,
SC_in_ThisMonth,
SC_in_ThisDay,
SC_in_GMToffset,
SC_out_StrBuffer,
SC_out_NoonTime_Hours,
SC_out_NoonTime_Minutes)
GetLocalSunSetTime (SC_in_Latitude,
SC_in_Longitude,
SC_in_ThisYear,
SC_in_ThisMonth,
SC_in_ThisDay,
SC_in_GMToffset,
SC_out_StrBuffer,
SC_out_SetTime_Hours,
SC_out_SetTime_Minutes)
' Free running cycle led
'-----------------------
if LED6 = false then
LED6 = true
else
LED6 = false
end if
wend
end.
