Merge pull request #295 from dickermoshe/meeus-noaa-rise-set-fallback

NOAACalculator, MeeusCalculator: don't miss zmanim the sun only just reaches
Thanks @dickermoshe

Author: KosherJava

src/main/java/com/kosherjava/zmanim/util/MeeusCalculator.java

@@ -297,6 +297,18 @@ private double getSunRiseSetUTC(LocalDate localDate, double latitude, double lon
 		double equationOfTime = noon.equationOfTime;
 		double solarDeclination = noon.solarDeclination;
 		double hourAngle = getSunHourAngle(latitude, solarDeclination, zenith, solarEvent);
+		if (Double.isNaN(hourAngle)) {
+			// Using the Sun's declination at noon, the Sun never reaches this zenith today. But the declination changes
+			// during the day (by up to about 0.2° in 12 hours), so the Sun may still just reach it around midnight, when it
+			// is at its lowest. Before reporting that there is no event, try again with the declination and equation of time
+			// at that solar midnight: the one before sunrise, or the one after sunset.
+			double noonmin = getSolarNoonMidnightUTC(julianDay, longitude, SolarEvent.NOON);
+			double midnight = noonmin + (solarEvent == SolarEvent.SUNRISE ? -720 : 720);
+			double tmidnight = getJulianCenturiesFromJulianDay(julianDay + midnight / 1440.0);
+			equationOfTime = getEquationOfTime(tmidnight);
+			solarDeclination = getSunDeclination(tmidnight);
+			hourAngle = getSunHourAngle(latitude, solarDeclination, zenith, solarEvent);
+		}
 		double delta = longitude - Math.toDegrees(hourAngle);
 		double timeUTC = 720 + (4 * delta) - equationOfTime;
 

src/main/java/com/kosherjava/zmanim/util/NOAACalculator.java

@@ -475,6 +475,17 @@ private static double getSunRiseSetUTC(LocalDate localDate, double latitude, dou
 		double equationOfTime = getEquationOfTime(tnoon);
 		double solarDeclination = getSunDeclination(tnoon);
 		double hourAngle = getSunHourAngle(latitude, solarDeclination, zenith, solarEvent);
+		if (Double.isNaN(hourAngle)) {
+			// Using the Sun's declination at noon, the Sun never reaches this zenith today. But the declination changes
+			// during the day (by up to about 0.2° in 12 hours), so the Sun may still just reach it around midnight, when it
+			// is at its lowest. Before reporting that there is no event, try again with the declination and equation of time
+			// at that solar midnight: the one before sunrise, or the one after sunset.
+			double midnight = noonmin + (solarEvent == SolarEvent.SUNRISE ? -720 : 720);
+			double tmidnight = getJulianCenturiesFromJulianDay(julianDay + midnight / 1440.0);
+			equationOfTime = getEquationOfTime(tmidnight);
+			solarDeclination = getSunDeclination(tmidnight);
+			hourAngle = getSunHourAngle(latitude, solarDeclination, zenith, solarEvent);
+		}
 		double delta = longitude - Math.toDegrees(hourAngle);
 		double timeDiff = 4 * delta;
 		double timeUTC = 720 + timeDiff - equationOfTime;

src/test/java/com/kosherjava/zmanim/util/MeeusCalculatorTest.java

@@ -244,6 +244,30 @@ public void getUTCSunriseWithDeltaTDisabled() {
 		}
 	}
 
+	/**
+	 * Days when the Sun only just reaches the zenith, around midnight. Using the Sun's declination at noon the zenith looks
+	 * out of reach, but using the declination at solar midnight it is reached, so these return a time rather than NaN.
+	 */
+	@Test
+	public void getUTCSunRiseSetWhenOnlyJustReached() {
+		GeoLocation london = new GeoLocation("London", 51.5074, -0.1278, 0, TestLocations.UTC);
+		GeoLocation arctic = new GeoLocation("Arctic", 75.51550698948259, -0.9961193842273985, 1516.0617032522935,
+				TestLocations.UTC);
+		Object[][] cases = {
+			{"London 16.1\u00b0 alos, Sun just reaches it", london, "2024-06-03", 106.1, true, 0.17268715},
+			{"London 18\u00b0 tzais, Sun just reaches it", london, "2024-07-20", 108.0, false, 23.85897993},
+			{"75.5\u00b0N sunrise just after midnight", arctic, "2093-04-22", ZENITH, true, 0.22098753},
+		};
+		for (Object[] c : cases) {
+			GeoLocation location = (GeoLocation) c[1];
+			double zenith = (double) c[3];
+			double actual = (boolean) c[4]
+					? calculator().getUTCSunrise(date((String) c[2]), location, zenith, zenith == ZENITH)
+					: calculator().getUTCSunset(date((String) c[2]), location, zenith, zenith == ZENITH);
+			assertUtcHour((String) c[0], (double) c[5], actual);
+		}
+	}
+
 	@Test(expected = IllegalArgumentException.class)
 	public void getTimeAtAzimuthRejectsUnsupportedAzimuth() {
 		calculator().getTimeAtAzimuth(date("2017-10-17"), geo("NJ"), 123.0);

src/test/java/com/kosherjava/zmanim/util/NOAACalculatorTest.java

@@ -230,6 +230,30 @@ public void getUTCSunriseWithoutElevationAdjustment() {
 		}
 	}
 
+	/**
+	 * Days when the Sun only just reaches the zenith, around midnight. Using the Sun's declination at noon the zenith looks
+	 * out of reach, but using the declination at solar midnight it is reached, so these return a time rather than NaN.
+	 */
+	@Test
+	public void getUTCSunRiseSetWhenOnlyJustReached() {
+		GeoLocation london = new GeoLocation("London", 51.5074, -0.1278, 0, TestLocations.UTC);
+		GeoLocation arctic = new GeoLocation("Arctic", 75.51550698948259, -0.9961193842273985, 1516.0617032522935,
+				TestLocations.UTC);
+		Object[][] cases = {
+			{"London 16.1\u00b0 alos, Sun just reaches it", london, "2024-06-03", 106.1, true, 0.17183556},
+			{"London 18\u00b0 tzais, Sun just reaches it", london, "2024-07-20", 108.0, false, 23.85934340},
+			{"75.5\u00b0N sunrise just after midnight", arctic, "2093-04-22", ZENITH, true, 0.21591651},
+		};
+		for (Object[] c : cases) {
+			GeoLocation location = (GeoLocation) c[1];
+			double zenith = (double) c[3];
+			double actual = (boolean) c[4]
+					? calculator().getUTCSunrise(date((String) c[2]), location, zenith, zenith == ZENITH)
+					: calculator().getUTCSunset(date((String) c[2]), location, zenith, zenith == ZENITH);
+			assertUtcHour((String) c[0], (double) c[5], actual);
+		}
+	}
+
 	@Test(expected = IllegalArgumentException.class)
 	public void getTimeAtAzimuthRejectsUnsupportedAzimuth() {
 		calculator().getTimeAtAzimuth(date("2017-10-17"), geo("NJ"), 123.0);