Hello.
This script gives you a report over the First day of the seasons of a year, together with the dates of the Equnioxes and Solstices.
You’ll need Satimage.osax, in order to run this.
Enjoy.
script solarEvents
property scriptTitle : "Solar Equinoxes, Solstices and Quarters"
# Compute Quarters:
property quarterNames : {"First day of Spring", "First Day of Summer", "First Day of Autumn", "First Day of Winter"}
property equinoxNames : {"Vernal Equinox", "Summer Solstice", "Auturnal Equinox", "Winter Solstice"}
property quarterAngles : {315.0, 45.0, 135.0, 225.0}
property debug : true
property tlvl : me
property EpochJD2000 : 2.451545E+6
property JulianDaysInCentury : 36525
script solaris
# Solar coordinates low accuracy Meeus chapter 24 p. 151
on longitude(JDE_t)
local T_, L0, m, ecc, equOfC, trueLongitude, trueAnamoly, R, omega, apparent, meanObliguityOfEcliptic
set T_ to JulianCenturies of tlvl by JDE_t
set L0 to mapWithinModulus(geometricMeanLongitudeOfSun(T_), 360)
set m to mapWithinModulus(meanAnamolyOfSun(T_), 360)
set ecc to eccOfEarthOrbit(T_)
set equOfC to sunsEquationOfCenter(T_, m, L0)
set trueLongitude to mapWithinModulus((L0 + equOfC), 360)
set trueAnamoly to mapWithinModulus((m + equOfC), 360)
set R to sunsRadiusVector(ecc, trueAnamoly)
# Meeus (24.5)
set omega to mapWithinModulus((125.04 - (T_ * 1934.136)), 360)
set apparentLong to (trueLongitude - 0.00569 - (0.00478 * (sin (omega * pi / 180))))
# The suns right ascension, and declination, in equatorial coordinates
set meanObliguityOfEcliptic to (meanObliquityOfEcliptic_Laskar for JDE_t)
return mapWithinModulus((sunsRightAscension(meanObliguityOfEcliptic, apparentLong) * 180 / pi), 360)
end longitude
on geometricMeanLongitudeOfSun(T)
# Meeus (24.2)
return (280.46645 + (T * (3.600076983E+4 + (T * 3.032E-4))))
end geometricMeanLongitudeOfSun
on meanAnamolyOfSun(T)
# Meeus ( 24.3)
return (357.52919 + (T * (3.5999053E+4 + (T * (-1.559E-4 + (T * -4.8E-7))))))
end meanAnamolyOfSun
on eccOfEarthOrbit(T)
# Meeus (24.4)
return (0.016708617 + (T * (-4.2037E-5 - (T * 1.236E-7))))
end eccOfEarthOrbit
on sunsEquationOfCenter(T, m, L0)
return (((1.9146 - (T * (0.004817 - (T * 1.4E-5)))) * (sin (m * pi / 180))) + ((0.01993 - (T * 1.01E-4)) * (sin (2 * (m * pi / 180)))) + ((2.9E-4 * (sin (3 * (m * pi / 180))))))
end sunsEquationOfCenter
on sunsRadiusVector(eccentricity, trueAnamoly)
# Meeus (24.5)
return ((1.00000010018 * (1 - (eccentricity ^ 2))) / (1 + (eccentricity * (cos (trueAnamoly * pi / 180)))))
end sunsRadiusVector
on sunsRightAscension(trueObliguity, longitude)
# (24.6)
return atan2 {((cos (trueObliguity * pi / 180)) * (sin (longitude * pi / 180))), (cos (longitude * pi / 180))}
end sunsRightAscension
on sunsDeclination(trueObliguity, longitude)
# (24.7)
return (asin ((sin (trueObliguity * pi / 180)) * (sin (longitude * pi / 180))))
end sunsDeclination
on meanObliquityOfEcliptic_Laskar for JD
# Meuus (21.3) p.135
# Extracted from sidereal Lib.
local U_
set U_ to (JD - 2451545) / 3652500
return (23.439291111111 ¬
+ (U_ * (-1.300258333333 + (U_ * (-4.305555555556E-4 + (U_ * (0.555347222222 + ¬
(U_ * (-0.014272222222 + (U_ * (0.069352777778 + (U_ * (-0.010847222222 + ¬
(U_ * (0.001977777778 + (U_ * (0.007741666667 + (U_ * (0.001608333333 + ¬
(6.805555555556E-4 * U_))))))))))))))))))))
end meanObliquityOfEcliptic_Laskar
on mapWithinModulus(aDecNumber, aModulus)
set aDecNumber to aDecNumber mod aModulus
if aDecNumber < 0 then set aDecNumber to aDecNumber + aModulus
return aDecNumber
end mapWithinModulus
end script
script equinoxlib
property s_table : {¬
{485, 324.96, 1934.136}, ¬
{203, 337.23, 3.2964467E+4}, ¬
{199, 342.08, 20.186}, ¬
{182, 27.85, 4.45267112E+5}, ¬
{156, 73.14, 4.5036886E+4}, ¬
{136, 171.52, 2.2518443E+4}, ¬
{77, 222.54, 6.5928934E+4}, ¬
{74, 296.72, 3034.906}, ¬
{70, 243.58, 9037.513}, ¬
{58, 119.81, 3.3718147E+4}, ¬
{52, 297.17, 150.678}, ¬
{50, 21.02, 2281.226}, ¬
{45, 247.54, 2.9929562E+4}, ¬
{44, 325.15, 3.1555956E+4}, ¬
{29, 60.93, 4443.417}, ¬
{18, 155.12, 6.7555328E+4}, ¬
{17, 288.79, 4562.452}, ¬
{16, 198.04, 6.2894029E+4}, ¬
{14, 199.76, 3.1436921E+4}, ¬
{12, 95.39, 1.4577848E+4}, ¬
{12, 287.11, 3.1931756E+4}, ¬
{12, 320.81, 3.4777259E+4}, ¬
{9, 227.73, 1222.114}, ¬
{8, 15.45, 1.6859074E+4}}
on sumOfS(T)
local sigma
set sigma to 0
repeat with i from 1 to (length of s_table)
set sigma to sigma + ((item 1 of item i of s_table) * (cos (((item 2 of item i of s_table) * pi / 180) + ((item 3 of item i of s_table) * pi / 180 * T))))
end repeat
return sigma
end sumOfS
on B_year(YR)
return ((YR - 2000) / 1000)
end B_year
on jdeforvernalequinox(YR)
# Returns JDE, that must be corrected to mean time, and civiltime.
# the correction for civiltime, is done by adding 0.5 to the hours.
local JDE, T_, w, deltaLambda, s
set JDE to ApproxVernalEquinox(YR)
set T_ to (JulianCenturies by JDE) of tlvl
set w to (3.5999373E+4 * T_ - 2.47) * pi / 180
set deltaLambda to 1 + (0.00334 * (cos (w))) + 7.0E-4 * (cos (2 * w))
set s to sumOfS(T_)
return (JDE + (1.0E-5 * s) / deltaLambda)
end jdeforvernalequinox
on ApproxVernalEquinox(YR)
local JDE0
if YR ≤ 1000 then
set YR to (YR / 1000)
set JDE0 to 1.72113929189E+6 + (YR * (3.652421374E+5 + (YR * (0.06134 + (YR * (0.00111 + (YR * 7.1E-4)))))))
else
set YR to ((YR - 2000) / 1000)
set JDE0 to 2.45162380984E+6 + (YR * (3.6524237404E+5 + (YR * (0.05169 + (YR * (-0.00411 + (YR * -5.7E-4)))))))
end if
return JDE0
end ApproxVernalEquinox
on JDEForSummerSolstice(YR)
# Returns JDE, that must be corrected to mean time, and civiltime.
# the correction for civiltime, is done by adding 0.5 to the hours.
local JDE, T_, w, deltaLambda, s
set JDE to ApproxSummerSolstice(YR)
set T_ to (JulianCenturies by JDE) of tlvl
set w to (3.5999373E+4 * T_ - 2.47) * pi / 180
set deltaLambda to 1 + (0.00334 * (cos (w))) + 7.0E-4 * (cos (2 * w))
set s to sumOfS(T_)
return (JDE + (1.0E-5 * s) / deltaLambda)
end JDEForSummerSolstice
on ApproxSummerSolstice(YR)
local JDE0
if YR ≤ 1000 then
set YR to (YR / 1000)
set JDE0 to 1.72123325401E+6 + (YR * (3.6524172562E+5 + (YR * (-0.05323 + (YR * (0.00907 + (YR * 2.5E-4)))))))
else
set YR to ((YR - 2000) / 1000)
set JDE0 to 2.45171656767E+6 + (YR * (3.6524162603E+5 + (YR * (0.00325 + (YR * (0.00888 + (YR * -3.0E-4)))))))
end if
return JDE0
end ApproxSummerSolstice
on JDEForAuturnalEquinox(YR)
# Returns JDE, that must be corrected to mean time, and civiltime.
# the correction for civiltime, is done by adding 0.5 to the hours.
local JDE, T_, w, deltaLambda, s
set JDE to ApproxAuturnalEquinox(YR)
set T_ to (JulianCenturies by JDE) of tlvl
set w to (3.5999373E+4 * T_ - 2.47) * pi / 180
set deltaLambda to 1 + (0.00334 * (cos (w))) + 7.0E-4 * (cos (2 * w))
set s to sumOfS(T_)
return (JDE + (1.0E-5 * s) / deltaLambda)
end JDEForAuturnalEquinox
on ApproxAuturnalEquinox(YR)
local JDE0
if YR ≤ 1000 then
set YR to (YR / 1000)
set JDE0 to 1.72132570455E+6 + (YR * (3.6524249558E+5 + (YR * (-0.11677 + (YR * (-0.00297 + (YR * 7.4E-4)))))))
else
set YR to ((YR - 2000) / 1000)
set JDE0 to 2.45181021715E+6 + (YR * (3.6524201767E+5 + (YR * (-0.11575 + (YR * (0.00337 + (YR * 7.8E-4)))))))
end if
return JDE0
end ApproxAuturnalEquinox
on jdeforwintersolstice(YR)
# Returns JDE, that must be corrected to mean time, and civiltime.
# the correction for civiltime, is done by adding 0.5 to the hours.
local JDE, T_, w, deltaLambda, s
set JDE to ApproxWinterSolstice(YR)
set T_ to (JulianCenturies by JDE) of tlvl
set w to (3.5999373E+4 * T_ - 2.47) * pi / 180
set deltaLambda to 1 + (0.00334 * (cos (w))) + 7.0E-4 * (cos (2 * w))
set s to sumOfS(T_)
return (JDE + (1.0E-5 * s) / deltaLambda)
end jdeforwintersolstice
on ApproxWinterSolstice(YR)
local JDE0
if YR ≤ 1000 then
set YR to (YR / 1000)
set JDE0 to 1.72141439987E+6 + (YR * (3.6524288257E+5 + (YR * (-0.00769 + (YR * (-0.00933 - (YR * 6.0E-5)))))))
else
set YR to ((YR - 2000) / 1000)
set JDE0 to 2.45190005952E+6 + (YR * 3.6524274049E+5 + (YR * (-0.06223 + (YR * (-0.00823 + (YR * 3.2E-4))))))
end if
return JDE0
end ApproxWinterSolstice
end script
on run
set suggestion to (current date)'s year as string
set wanted_year to inputYear(suggestion)
set {equinoxTable, quarterTable} to {{}, {}}
# We do compute the previous Winter Solstice.
# We also keep all times in UT, and converts at the end.
set end of equinoxTable to my equinoxlib's jdeforwintersolstice(wanted_year - 1)
set end of equinoxTable to my equinoxlib's jdeforvernalequinox(wanted_year)
set end of equinoxTable to my equinoxlib's JDEForSummerSolstice(wanted_year)
set end of equinoxTable to my equinoxlib's JDEForAuturnalEquinox(wanted_year)
set end of equinoxTable to my equinoxlib's jdeforwintersolstice(wanted_year)
repeat with q from 1 to 4
# compute Arithmetic Mean and prepare delta
set JDMean to ((item q of equinoxTable) + (item (q + 1) of equinoxTable)) / 2
set {delta, n} to {60, 1}
# Use binary search, to iterate towards the solution.
repeat while true
set curLong to my solaris's longitude(JDMean)
if (abs ((item q of quarterAngles) - curLong)) < 0.1 then
set end of quarterTable to JDMean
exit repeat
end if
if curLong > (item q of quarterAngles) then
set JDMean to JDMean - (delta / (n * 2))
else
set JDMean to JDMean + (delta / (n * 2))
end if
set n to n + 1
end repeat
end repeat
_MAKEDOC((scriptTitle & " For " & wanted_year as text) & "." & linefeed & linefeed, 600, 300)
repeat with i from 1 to 8
local yr_, mnt_, day_, h_, m_, s_
if (i mod 2 = 1) then # odd - Quarter!
set {yr_, mnt_, day_, h_, m_, s_} to adjustJulianDateAccordingToGMT(item (i div 2 + 1) of quarterTable)
_PRINT(formatDate(yr_, mnt_, day_) & " : " & item (i div 2 + 1) of quarterNames & linefeed)
else
set {yr_, mnt_, day_, h_, m_, s_} to adjustJulianDateAccordingToGMT(item (i div 2 + 1) of equinoxTable)
_PRINT(formatDate(yr_, mnt_, day_) & " : " & item (i div 2) of equinoxNames & linefeed)
end if
end repeat
end run
on TZtoTZ(TZ1date, TZ1, TZ2)
return (do shell script ("eraTime=$(TZ=" & TZ1 & " date -jf '%FT%T' '" & (TZ1date as «class isot» as string) & "' '+%s') ; TZ=" & TZ2 & " date -r \"$eraTime\" '+%FT%T'") as «class isot») as date
end TZtoTZ
on civilTimeToGmtFromLocalTZInSeconds(anASDate)
return (anASDate - TZtoTZ(anASDate, (do shell script "readlink /etc/localtime |sed -n 's_/[^/]*/[^/]*/[^/]*/\\([^/]*\\)/\\([^/]*\\).*_\\1/\\2_p'"), "GMT"))
end civilTimeToGmtFromLocalTZInSeconds
on lastdayInMonth(YR, mo)
local daynumbers
set daynumbers to {31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31}
if mo < 1 or mo > 12 then
error "lastDay: Month out of range..." number 3010
else if mo = 2 then
if ((YR mod 4) = 0) and (((YR mod 100) ≠0) or ((YR mod 400) = 0)) and (mo > 2) then
set item 2 of daynumbers to 29
end if
end if
return item mo of daynumbers
end lastdayInMonth
on set_date_long(YR, mnth, dy, hr, min, sec)
# NG http://macscripter.net/viewtopic.php?id=41374 there is a short version there too.
tell (current date) to set {day, year, its month, day, its hours, its minutes, its seconds, cur_date} to {1, YR, mnth, dy, hr, min, sec, it}
return cur_date
end set_date_long
# returns a time triplet, from a Julian day number,
# Adjusted to return UT.
on timeTriplet from JDN
# This returns as standard time triplet, where the first point of time of a day is 00:00
# and the last point of a time is 23:59:59.9999999 and so on...
local timeFraction, tempvar, _hours, _minutes, _seconds
set timeFraction to (JDN mod 1) * 24
set tempvar to (timeFraction div 1 as integer)
set _hours to (tempvar + 12) mod 24
set timeFraction to (timeFraction - tempvar)
set _minutes to timeFraction * 60 div 1
set _seconds to (timeFraction - (_minutes / 60)) * 3600
return {_hours, _minutes, _seconds}
end timeTriplet
on adjustJulianDateAccordingToGMT(JDE)
local dt, _yr, mnt, _day, _hr, _min, _sec
# returns a sextet
set JDE to JDE + 0.5
# Find day and time of the event so far in ut time.
set {_yr, mnt, _day} to gregorianCalendarTriplet from JDE
set {_hr, _min, _sec} to timeTriplet from JDE
# we have to adjust the hours
if _hr > 12 then
set _hr to _hr - 12
else
set _hr to _hr + 12
end if
# figure out the TimeToGMT at this instace
set dt to set_date_long(_yr, mnt, _day, _hr, _min, _sec)
set _hr to _hr + (civilTimeToGmtFromLocalTZInSeconds(dt) / 3600)
if _hr > 24 then
set _hr to _hr - 24
set _day to _day + 1
if not _day ≤ lastdayInMonth(_yr, mnt) then
if mnt < 12 then
set _day to 1
set mnt to mnt + 1
else
set _yr to _yr + 1
set mnt to 1
set _day to 1
end if
end if
else if _hr < 0 then
set _hr to _hr + 24
if _day > 1 then
set _day to _day - 1
else if mnt > 1 then
set mnt to mnt - 1
set _day to lastdayInMonth(_yr, mnt)
else
set _yr to _yr - 1
set mnt to 12
set _day to 31
end if
end if
return {_yr, mnt, _day, _hr, _min, _sec}
end adjustJulianDateAccordingToGMT
# returns the number of Julian centuries that has passed
# since the JDE
on JulianCenturies by JDE
# Meeus (11.1)
return ((JDE - (my EpochJD2000)) / (my JulianDaysInCentury))
end JulianCenturies
on gregorianCalendarTriplet from JDN
# This is an algorithm by Richards to convert a Julian Day Number, J, to a date in the Gregorian calendar
# (proleptic, when applicable). Richards does not state which dates the algorithm is valid for.
# it seems to work for dates after, BUT NOT INCLUDING 15th of oct 1582.
# I have extended it to return a hour, min, sec triplet for the fractional part of the JDN.
local y, m, n, R, p, V, s, w, b, c, f, g, h, D, MT, YR, JD
set y to 4716
set j to 1401
set m to 2
set n to 12
set R to 4
set p to 1461
set V to 3
set u to 5
set s to 153
set w to 2
set b to 274277
set c to -38
set JD to JDN div 1
set f to JD + j + (((4 * JD + b) div 146097) * 3) div 4 + c
set e to R * f + V
set g to (e mod p) div R
set h to u * g + w
set D to (h mod s) div u + 1
set MT to (h div s + m) mod n + 1
set YR to e div p - y + (n + m - MT) div n
return {YR, MT, D}
end gregorianCalendarTriplet
on _PRINT(txt)
# Nigel Garvey made the original version.
local a
tell application "TextEdit"
if not my debug then
activate
else
do shell script "open -b \"com.apple.textedit\""
end if
if ((count its documents) is 0) then make new document
set a to (get path of its front document)
try
if a is "" then set a to missing value
on error
set a to missing value
end try
if a is not missing value then make new document
make new paragraph at end of text of front document with data (txt & linefeed) with properties {font:"Andale Mono", size:12}
end tell
return true
end _PRINT
on formatDate(YR, mnt, _day)
global us_citizen
local a
try
set a to us_citizen
on error
set us_citizen to getUsCitizenShip()
end try
if us_citizen then
return (pad(mnt) & "/" & pad(_day) & "/" & YR)
else
return (pad(_day) & "/" & pad(mnt) & "/" & YR)
end if
end formatDate
on pad(n)
-- 22. September 2013, improved pad, for the cases where n < 0 and n > 100.
-- Thanks to Yvan Koenig.
-- Fixed a bug
if n < -100 then
return n as text
else if n < 0 then
return ("-" & text -2 thru -1 of ("0" & ((n as number) * -1) as text))
-- error "The pad handler assumes positive values" number 5000
else if n < 100 then
return text 2 thru -1 of ((100 + n) as text)
else
return n as text
end if
end pad
to getUsCitizenShip()
try
set localeString to (do shell script "defaults read .GlobalPreferences AppleLocale 2>/dev/null || echo en_US |tr '
' ' '")
on error e number n
return true
end try
if localeString is "en_US" then return true
return false
end getUsCitizenShip
on _MAKEDOC(heading, _width, _height)
tell application "TextEdit"
tell (make new document at front)
set text of it to heading
set font of every character to "Andale Mono"
set size of every character to 13
end tell
tell its front window
# set {item 3 of its bounds, item 4 of its bounds} to {800, 300}
set _bounds to bounds
set item 3 of _bounds to (item 1 of _bounds) + _width
set item 4 of _bounds to (item 2 of _bounds) + _height
set bounds to _bounds
end tell
end tell
end _MAKEDOC
on inputYear(suggestion)
local stdMsg, probe, theMessage, anicon
set stdMsg to "Enter the year you want to see the Solar Phases for."
set theMessage to stdMsg
set anicon to note
repeat while true
tell application (path to frontmost application as text)
try
set theDigitstring to text returned of (display dialog theMessage default answer suggestion with title my scriptTitle with icon anicon buttons {"Quit", "Ok"} cancel button 1 default button 2)
on error
error number -128
end try
end tell
try
set probe to theDigitstring as number
exit repeat
on error
set theMessage to "\"" & theDigit & "\" isn't a valid year!" & return & stdMsg
set anicon to caution
end try
end repeat
return theDigitstring
end inputYear
end script
tell solarEvents to run