From 72d154c5e09a997fc6cb16700594fa28d2bea971 Mon Sep 17 00:00:00 2001
From: Dianne Skoll <dfs@roaringpenguin.com>
Date: Thu, 19 Jan 2017 10:11:54 -0500
Subject: [PATCH] Patch from Stephen Morgan to calculate astronomical and  
 nautical twilight in addition to civil twilight.

---
 docs/WHATSNEW |  5 +++++
 man/remind.1  | 18 ++++++++++++++++++
 src/funcs.c   | 46 +++++++++++++++++++++++++++++++++++++++-------
 3 files changed, 62 insertions(+), 7 deletions(-)

diff --git a/docs/WHATSNEW b/docs/WHATSNEW
index 0ba78f43..3e00dce3 100644
--- a/docs/WHATSNEW
+++ b/docs/WHATSNEW
@@ -1,5 +1,10 @@
 CHANGES TO REMIND
 
+* Version 3.1 Patch 16 - 2017-XX-XX
+
+- IMPROVEMENT: Add patch from Stephen Morgan to calculate astronomical and
+  nautical twilight in addition to civil twilight.
+
 * Version 3.1 Patch 15 - 2015-07-27
 
 - BUG FIX: Fix a buffer overflow found by Alexander Keller
diff --git a/man/remind.1 b/man/remind.1
index e430dfe4..74425ffa 100644
--- a/man/remind.1
+++ b/man/remind.1
@@ -2123,6 +2123,15 @@ be a number as described in the UNIX \fBaccess\fR(2) system call.  The
 function returns 0 if the file can be accessed with the specified \fImode\fR,
 and \-1 otherwise.
 .TP
+.B adawn([dq_date])
+Returns the time of "astronomical dawn" on the specified \fIdate\fR.  If
+\fIdate\fR is omitted, defaults to \fBtoday()\fR.  If a \fIdatetime\fR object
+is supplied, only the date component is used.
+.TP
+.B adusk([dq_date])
+Returns the time of "astronomical twilight" on the specified \fIdate\fR.  If
+\fIdate\fR is omitted, defaults to \fBtoday()\fR.
+.TP
 .B args(s_fname)
 Returns the number of arguments expected by the user-defined function
 \fIfname\fR, or \-1 if no such user-defined function exists.  Note that
@@ -2517,6 +2526,15 @@ which default to \fBtoday()\fR and midnight, respectively.  The returned
 value is an integer from 0 to 359, representing the phase of the moon
 in degrees.  0 is a new moon, 180 is a full moon, 90 is first-quarter, etc.
 .TP
+.B ndawn([dq_date])
+Returns the time of "nautical dawn" on the specified \fIdate\fR.  If
+\fIdate\fR is omitted, defaults to \fBtoday()\fR.  If a \fIdatetime\fR object
+is supplied, only the date component is used.
+.TP
+.B ndusk([dq_date])
+Returns the time of "nautical twilight" on the specified \fIdate\fR.  If
+\fIdate\fR is omitted, defaults to \fBtoday()\fR.
+.TP
 .B nonomitted(dq_start, dq_end [,s_wkday...])
 This function returns the number of \fInon-\fRomitted days between
 \fIstart\fR and \fIend\fR.  If \fIstart\fR is non-omitted, then it is
diff --git a/src/funcs.c b/src/funcs.c
index e8347891..fbf60662 100644
--- a/src/funcs.c
+++ b/src/funcs.c
@@ -110,6 +110,10 @@ static int	FPsshade	(func_info *);
 static	int	FShell		(func_info *);
 static	int	FStrlen		(func_info *);
 static	int	FSubstr		(func_info *);
+static	int	FADawn		(func_info *);
+static	int	FADusk	 	(func_info *);
+static	int	FNDawn		(func_info *);
+static	int	FNDusk	 	(func_info *);
 static	int	FDawn		(func_info *);
 static	int	FDusk	 	(func_info *);
 static	int	FSunset		(func_info *);
@@ -188,6 +192,8 @@ BuiltinFunc Func[] = {
 
     {   "abs",          1,      1,      1,          FAbs },
     {   "access",       2,      2,      0,          FAccess },
+    {   "adawn",        0,      1,      0,          FADawn},
+    {   "adusk",        0,      1,      0,          FADusk},
     {   "args",         1,      1,      0,          FArgs   },
     {   "asc",          1,      1,      1,          FAsc    },
     {   "baseyr",       0,      0,      1,          FBaseyr },
@@ -233,6 +239,8 @@ BuiltinFunc Func[] = {
     {   "moondatetime", 1,      3,      0,          FMoondatetime },
     {   "moonphase",    0,      2,      0,          FMoonphase },
     {   "moontime",     1,      3,      0,          FMoontime },
+    {   "ndawn",        0,      1,      0,          FNDawn},
+    {   "ndusk",        0,      1,      0,          FNDusk},
     {   "nonomitted",   2,      NO_MAX, 0,          FNonomitted },
     {   "now",          0,      0,      0,          FNow    },
     {   "ord",          1,      1,      1,          FOrd    },
@@ -1878,9 +1886,6 @@ static int SunStuff(int rise, double cosz, int jul)
 
     FromJulian(jul, &year, &mon, &day);
 
-    if (rise > 1)
-	rise -= 2;
-
 /* Following formula on page B6 exactly... */
     t = (double) jul;
     if (rise) {
@@ -1970,19 +1975,28 @@ static int SunStuff(int rise, double cosz, int jul)
 static int FSun(int rise, func_info *info)
 {
     int jul = JulianToday;
-    double cosz = -0.014543897;  /* for sunrise and sunset */
+    double cosz;
     int r;
 
+    if (rise == 0 || rise == 1) {
+    /* Sunrise and sunset : cos(90 degrees + 50 arcminutes) */
+    cosz = -0.01454389765158243;
+    } else if (rise == 2 || rise == 3) {
     /* Civil twilight: cos(96 degrees) */
-    if (rise == 2 || rise == 3) {
-	cosz = -0.104528463268;
+	cosz = -0.10452846326765333;
+    } else if (rise == 4 || rise == 5) {
+    /* Nautical twilight: cos(102 degrees) */
+	cosz = -0.20791169081775912;
+    } else if (rise == 6 || rise == 7) {
+    /* Astronomical twilight: cos(108 degrees) */
+	cosz = -0.30901699437494734;
     }
     if (Nargs >= 1) {
 	if (!HASDATE(ARG(0))) return E_BAD_TYPE;
 	jul = DATEPART(ARG(0));
     }
 
-    r = SunStuff(rise, cosz, jul);
+    r = SunStuff(rise % 2, cosz, jul);
     if (r == NO_TIME) {
 	RETVAL = 0;
 	RetVal.type = INT_TYPE;
@@ -2014,6 +2028,24 @@ static int FDusk(func_info *info)
     return FSun(2, info);
 }
 
+static int FNDawn(func_info *info)
+{
+    return FSun(5, info);
+}
+static int FNDusk(func_info *info)
+{
+    return FSun(4, info);
+}
+
+static int FADawn(func_info *info)
+{
+    return FSun(7, info);
+}
+static int FADusk(func_info *info)
+{
+    return FSun(6, info);
+}
+
 /***************************************************************/
 /*                                                             */
 /*  FFiledate                                                  */