| | | 1 | | // Licensed to the .NET Foundation under one or more agreements. |
| | | 2 | | // The .NET Foundation licenses this file to you under the MIT license. |
| | | 3 | | |
| | | 4 | | using System.Diagnostics; |
| | | 5 | | |
| | | 6 | | namespace System.Globalization |
| | | 7 | | { |
| | | 8 | | internal static class CalendricalCalculationsHelper |
| | | 9 | | { |
| | | 10 | | private const double FullCircleOfArc = 360.0; // 360.0; |
| | | 11 | | private const int HalfCircleOfArc = 180; |
| | | 12 | | private const double TwelveHours = 0.5; // half a day |
| | | 13 | | private const double Noon2000Jan01 = 730120.5; |
| | | 14 | | internal const double MeanTropicalYearInDays = 365.242189; |
| | | 15 | | private const double MeanSpeedOfSun = MeanTropicalYearInDays / FullCircleOfArc; |
| | | 16 | | private const double LongitudeSpring = 0.0; |
| | | 17 | | private const double TwoDegreesAfterSpring = 2.0; |
| | | 18 | | private const int DaysInUniformLengthCentury = 36525; |
| | | 19 | | |
| | 0 | 20 | | private static readonly long s_startOf1810 = GetNumberOfDays(new DateTime(1810, 1, 1)); |
| | 0 | 21 | | private static readonly long s_startOf1900Century = GetNumberOfDays(new DateTime(1900, 1, 1)); |
| | | 22 | | |
| | 0 | 23 | | private static ReadOnlySpan<double> Coefficients1900to1987 => [-0.00002, 0.000297, 0.025184, -0.181133, 0.553040 |
| | 0 | 24 | | private static ReadOnlySpan<double> Coefficients1800to1899 => [-0.000009, 0.003844, 0.083563, 0.865736, 4.867575 |
| | 0 | 25 | | private static ReadOnlySpan<double> Coefficients1700to1799 => [8.118780842, -0.005092142, 0.003336121, -0.000026 |
| | 0 | 26 | | private static ReadOnlySpan<double> Coefficients1620to1699 => [196.58333, -4.0675, 0.0219167]; |
| | 0 | 27 | | private static ReadOnlySpan<double> LambdaCoefficients => [280.46645, 36000.76983, 0.0003032]; |
| | 0 | 28 | | private static ReadOnlySpan<double> AnomalyCoefficients => [357.52910, 35999.05030, -0.0001559, -0.00000048]; |
| | 0 | 29 | | private static ReadOnlySpan<double> EccentricityCoefficients => [0.016708617, -0.000042037, -0.0000001236]; |
| | 0 | 30 | | private static ReadOnlySpan<double> CoefficientsA => [124.90, -1934.134, 0.002063]; |
| | 0 | 31 | | private static ReadOnlySpan<double> CoefficientsB => [201.11, 72001.5377, 0.00057]; |
| | 0 | 32 | | private static ReadOnlySpan<double> Coefficients => [23.43929111111111, -0.013004166666666667, -1.63888888888888 |
| | | 33 | | |
| | | 34 | | private static double RadiansFromDegrees(double degree) |
| | | 35 | | { |
| | 0 | 36 | | return degree * Math.PI / 180; |
| | | 37 | | } |
| | | 38 | | |
| | | 39 | | private static double SinOfDegree(double degree) |
| | | 40 | | { |
| | 0 | 41 | | return Math.Sin(RadiansFromDegrees(degree)); |
| | | 42 | | } |
| | | 43 | | |
| | | 44 | | private static double CosOfDegree(double degree) |
| | | 45 | | { |
| | 0 | 46 | | return Math.Cos(RadiansFromDegrees(degree)); |
| | | 47 | | } |
| | | 48 | | |
| | | 49 | | private static double TanOfDegree(double degree) |
| | | 50 | | { |
| | 0 | 51 | | return Math.Tan(RadiansFromDegrees(degree)); |
| | | 52 | | } |
| | | 53 | | |
| | | 54 | | private static double Obliquity(double julianCenturies) |
| | | 55 | | { |
| | 0 | 56 | | return PolynomialSum(Coefficients, julianCenturies); |
| | | 57 | | } |
| | | 58 | | |
| | | 59 | | internal static long GetNumberOfDays(DateTime date) |
| | | 60 | | { |
| | 0 | 61 | | return date.Ticks / TimeSpan.TicksPerDay; |
| | | 62 | | } |
| | | 63 | | |
| | | 64 | | private static int GetGregorianYear(double numberOfDays) |
| | | 65 | | { |
| | 0 | 66 | | return new DateTime(Math.Min((long)(Math.Floor(numberOfDays) * TimeSpan.TicksPerDay), DateTime.MaxValue.Tick |
| | | 67 | | } |
| | | 68 | | |
| | | 69 | | private enum CorrectionAlgorithm |
| | | 70 | | { |
| | | 71 | | Default, |
| | | 72 | | Year1988to2019, |
| | | 73 | | Year1900to1987, |
| | | 74 | | Year1800to1899, |
| | | 75 | | Year1700to1799, |
| | | 76 | | Year1620to1699 |
| | | 77 | | } |
| | | 78 | | |
| | | 79 | | private readonly struct EphemerisCorrectionAlgorithmMap |
| | | 80 | | { |
| | | 81 | | public EphemerisCorrectionAlgorithmMap(int year, CorrectionAlgorithm algorithm) |
| | | 82 | | { |
| | 0 | 83 | | _lowestYear = year; |
| | 0 | 84 | | _algorithm = algorithm; |
| | 0 | 85 | | } |
| | | 86 | | |
| | | 87 | | internal readonly int _lowestYear; |
| | | 88 | | internal readonly CorrectionAlgorithm _algorithm; |
| | | 89 | | } |
| | | 90 | | |
| | 0 | 91 | | private static readonly EphemerisCorrectionAlgorithmMap[] s_ephemerisCorrectionTable = |
| | 0 | 92 | | [ |
| | 0 | 93 | | // lowest year that starts algorithm, algorithm to use |
| | 0 | 94 | | new EphemerisCorrectionAlgorithmMap(2020, CorrectionAlgorithm.Default), |
| | 0 | 95 | | new EphemerisCorrectionAlgorithmMap(1988, CorrectionAlgorithm.Year1988to2019), |
| | 0 | 96 | | new EphemerisCorrectionAlgorithmMap(1900, CorrectionAlgorithm.Year1900to1987), |
| | 0 | 97 | | new EphemerisCorrectionAlgorithmMap(1800, CorrectionAlgorithm.Year1800to1899), |
| | 0 | 98 | | new EphemerisCorrectionAlgorithmMap(1700, CorrectionAlgorithm.Year1700to1799), |
| | 0 | 99 | | new EphemerisCorrectionAlgorithmMap(1620, CorrectionAlgorithm.Year1620to1699), |
| | 0 | 100 | | new EphemerisCorrectionAlgorithmMap(int.MinValue, CorrectionAlgorithm.Default) // default must be last |
| | 0 | 101 | | ]; |
| | | 102 | | |
| | | 103 | | private static double Reminder(double divisor, double dividend) |
| | | 104 | | { |
| | 0 | 105 | | double whole = Math.Floor(divisor / dividend); |
| | 0 | 106 | | return divisor - (dividend * whole); |
| | | 107 | | } |
| | | 108 | | |
| | | 109 | | private static double NormalizeLongitude(double longitude) |
| | | 110 | | { |
| | 0 | 111 | | longitude = Reminder(longitude, FullCircleOfArc); |
| | 0 | 112 | | if (longitude < 0) |
| | | 113 | | { |
| | 0 | 114 | | longitude += FullCircleOfArc; |
| | | 115 | | } |
| | 0 | 116 | | return longitude; |
| | | 117 | | } |
| | | 118 | | |
| | | 119 | | public static double AsDayFraction(double longitude) |
| | | 120 | | { |
| | 0 | 121 | | return longitude / FullCircleOfArc; |
| | | 122 | | } |
| | | 123 | | |
| | | 124 | | private static double PolynomialSum(ReadOnlySpan<double> coefficients, double indeterminate) |
| | | 125 | | { |
| | 0 | 126 | | double sum = coefficients[0]; |
| | 0 | 127 | | double indeterminateRaised = 1; |
| | 0 | 128 | | for (int i = 1; i < coefficients.Length; i++) |
| | | 129 | | { |
| | 0 | 130 | | indeterminateRaised *= indeterminate; |
| | 0 | 131 | | sum += (coefficients[i] * indeterminateRaised); |
| | | 132 | | } |
| | | 133 | | |
| | 0 | 134 | | return sum; |
| | | 135 | | } |
| | | 136 | | |
| | | 137 | | private static double CenturiesFrom1900(int gregorianYear) |
| | | 138 | | { |
| | 0 | 139 | | long july1stOfYear = GetNumberOfDays(new DateTime(gregorianYear, 7, 1)); |
| | 0 | 140 | | return (double)(july1stOfYear - s_startOf1900Century) / DaysInUniformLengthCentury; |
| | | 141 | | } |
| | | 142 | | |
| | | 143 | | // the following formulas defines a polynomial function which gives us the amount that the earth is slowing down |
| | | 144 | | private static double DefaultEphemerisCorrection(int gregorianYear) |
| | | 145 | | { |
| | 0 | 146 | | Debug.Assert(gregorianYear < 1620 || 2020 <= gregorianYear); |
| | 0 | 147 | | long january1stOfYear = GetNumberOfDays(new DateTime(gregorianYear, 1, 1)); |
| | 0 | 148 | | double daysSinceStartOf1810 = january1stOfYear - s_startOf1810; |
| | 0 | 149 | | double x = TwelveHours + daysSinceStartOf1810; |
| | 0 | 150 | | return ((Math.Pow(x, 2) / 41048480) - 15) / TimeSpan.SecondsPerDay; |
| | | 151 | | } |
| | | 152 | | |
| | | 153 | | private static double EphemerisCorrection1988to2019(int gregorianYear) |
| | | 154 | | { |
| | 0 | 155 | | Debug.Assert(1988 <= gregorianYear && gregorianYear <= 2019); |
| | 0 | 156 | | return (double)(gregorianYear - 1933) / TimeSpan.SecondsPerDay; |
| | | 157 | | } |
| | | 158 | | |
| | | 159 | | private static double EphemerisCorrection1900to1987(int gregorianYear) |
| | | 160 | | { |
| | 0 | 161 | | Debug.Assert(1900 <= gregorianYear && gregorianYear <= 1987); |
| | 0 | 162 | | double centuriesFrom1900 = CenturiesFrom1900(gregorianYear); |
| | 0 | 163 | | return PolynomialSum(Coefficients1900to1987, centuriesFrom1900); |
| | | 164 | | } |
| | | 165 | | |
| | | 166 | | private static double EphemerisCorrection1800to1899(int gregorianYear) |
| | | 167 | | { |
| | 0 | 168 | | Debug.Assert(1800 <= gregorianYear && gregorianYear <= 1899); |
| | 0 | 169 | | double centuriesFrom1900 = CenturiesFrom1900(gregorianYear); |
| | 0 | 170 | | return PolynomialSum(Coefficients1800to1899, centuriesFrom1900); |
| | | 171 | | } |
| | | 172 | | |
| | | 173 | | private static double EphemerisCorrection1700to1799(int gregorianYear) |
| | | 174 | | { |
| | 0 | 175 | | Debug.Assert(1700 <= gregorianYear && gregorianYear <= 1799); |
| | 0 | 176 | | double yearsSince1700 = gregorianYear - 1700; |
| | 0 | 177 | | return PolynomialSum(Coefficients1700to1799, yearsSince1700) / TimeSpan.SecondsPerDay; |
| | | 178 | | } |
| | | 179 | | |
| | | 180 | | private static double EphemerisCorrection1620to1699(int gregorianYear) |
| | | 181 | | { |
| | 0 | 182 | | Debug.Assert(1620 <= gregorianYear && gregorianYear <= 1699); |
| | 0 | 183 | | double yearsSince1600 = gregorianYear - 1600; |
| | 0 | 184 | | return PolynomialSum(Coefficients1620to1699, yearsSince1600) / TimeSpan.SecondsPerDay; |
| | | 185 | | } |
| | | 186 | | |
| | | 187 | | // ephemeris-correction: correction to account for the slowing down of the rotation of the earth |
| | | 188 | | private static double EphemerisCorrection(double time) |
| | | 189 | | { |
| | 0 | 190 | | int year = GetGregorianYear(time); |
| | 0 | 191 | | foreach (EphemerisCorrectionAlgorithmMap map in s_ephemerisCorrectionTable) |
| | | 192 | | { |
| | 0 | 193 | | if (map._lowestYear <= year) |
| | | 194 | | { |
| | 0 | 195 | | switch (map._algorithm) |
| | | 196 | | { |
| | 0 | 197 | | case CorrectionAlgorithm.Default: return DefaultEphemerisCorrection(year); |
| | 0 | 198 | | case CorrectionAlgorithm.Year1988to2019: return EphemerisCorrection1988to2019(year); |
| | 0 | 199 | | case CorrectionAlgorithm.Year1900to1987: return EphemerisCorrection1900to1987(year); |
| | 0 | 200 | | case CorrectionAlgorithm.Year1800to1899: return EphemerisCorrection1800to1899(year); |
| | 0 | 201 | | case CorrectionAlgorithm.Year1700to1799: return EphemerisCorrection1700to1799(year); |
| | 0 | 202 | | case CorrectionAlgorithm.Year1620to1699: return EphemerisCorrection1620to1699(year); |
| | | 203 | | } |
| | | 204 | | |
| | | 205 | | break; // break the loop and assert eventually |
| | | 206 | | } |
| | | 207 | | } |
| | | 208 | | |
| | 0 | 209 | | Debug.Fail("Not expected to come here"); |
| | | 210 | | return DefaultEphemerisCorrection(year); |
| | | 211 | | } |
| | | 212 | | |
| | | 213 | | public static double JulianCenturies(double moment) |
| | | 214 | | { |
| | 0 | 215 | | double dynamicalMoment = moment + EphemerisCorrection(moment); |
| | 0 | 216 | | return (dynamicalMoment - Noon2000Jan01) / DaysInUniformLengthCentury; |
| | | 217 | | } |
| | | 218 | | |
| | | 219 | | // equation-of-time; approximate the difference between apparent solar time and mean solar time |
| | | 220 | | // formal definition is EOT = GHA - GMHA |
| | | 221 | | // GHA is the Greenwich Hour Angle of the apparent (actual) Sun |
| | | 222 | | // GMHA is the Greenwich Mean Hour Angle of the mean (fictitious) Sun |
| | | 223 | | // http://www.esrl.noaa.gov/gmd/grad/solcalc/ |
| | | 224 | | // http://en.wikipedia.org/wiki/Equation_of_time |
| | | 225 | | private static double EquationOfTime(double time) |
| | | 226 | | { |
| | 0 | 227 | | double julianCenturies = JulianCenturies(time); |
| | 0 | 228 | | double lambda = PolynomialSum(LambdaCoefficients, julianCenturies); |
| | 0 | 229 | | double anomaly = PolynomialSum(AnomalyCoefficients, julianCenturies); |
| | 0 | 230 | | double eccentricity = PolynomialSum(EccentricityCoefficients, julianCenturies); |
| | | 231 | | |
| | 0 | 232 | | double epsilon = Obliquity(julianCenturies); |
| | 0 | 233 | | double tanHalfEpsilon = TanOfDegree(epsilon / 2); |
| | 0 | 234 | | double y = tanHalfEpsilon * tanHalfEpsilon; |
| | | 235 | | |
| | 0 | 236 | | double dividend = ((y * SinOfDegree(2 * lambda)) |
| | 0 | 237 | | - (2 * eccentricity * SinOfDegree(anomaly)) |
| | 0 | 238 | | + (4 * eccentricity * y * SinOfDegree(anomaly) * CosOfDegree(2 * lambda)) |
| | 0 | 239 | | - (0.5 * Math.Pow(y, 2) * SinOfDegree(4 * lambda)) |
| | 0 | 240 | | - (1.25 * Math.Pow(eccentricity, 2) * SinOfDegree(2 * anomaly))); |
| | | 241 | | const double Divisor = 2 * Math.PI; |
| | 0 | 242 | | double equation = dividend / Divisor; |
| | | 243 | | |
| | | 244 | | // approximation of equation of time is not valid for dates that are many millennia in the past or future |
| | | 245 | | // thus limited to a half day |
| | 0 | 246 | | return Math.CopySign(Math.Min(Math.Abs(equation), TwelveHours), equation); |
| | | 247 | | } |
| | | 248 | | |
| | | 249 | | private static double AsLocalTime(double apparentMidday, double longitude) |
| | | 250 | | { |
| | | 251 | | // slightly inaccurate since equation of time takes mean time not apparent time as its argument, but the dif |
| | 0 | 252 | | double universalTime = apparentMidday - AsDayFraction(longitude); |
| | 0 | 253 | | return apparentMidday - EquationOfTime(universalTime); |
| | | 254 | | } |
| | | 255 | | |
| | | 256 | | // midday |
| | | 257 | | public static double Midday(double date, double longitude) |
| | | 258 | | { |
| | 0 | 259 | | return AsLocalTime(date + TwelveHours, longitude) - AsDayFraction(longitude); |
| | | 260 | | } |
| | | 261 | | |
| | | 262 | | private static double InitLongitude(double longitude) |
| | | 263 | | { |
| | 0 | 264 | | return NormalizeLongitude(longitude + HalfCircleOfArc) - HalfCircleOfArc; |
| | | 265 | | } |
| | | 266 | | |
| | | 267 | | // midday-in-tehran |
| | | 268 | | public static double MiddayAtPersianObservationSite(double date) |
| | | 269 | | { |
| | 0 | 270 | | return Midday(date, InitLongitude(52.5)); // 52.5 degrees east - longitude of UTC+3:30 which defines Iranian |
| | | 271 | | } |
| | | 272 | | |
| | | 273 | | private static double PeriodicTerm(double julianCenturies, int x, double y, double z) |
| | | 274 | | { |
| | 0 | 275 | | return x * SinOfDegree(y + z * julianCenturies); |
| | | 276 | | } |
| | | 277 | | |
| | | 278 | | private static double SumLongSequenceOfPeriodicTerms(double julianCenturies) |
| | | 279 | | { |
| | 0 | 280 | | double sum = 0.0; |
| | 0 | 281 | | sum += PeriodicTerm(julianCenturies, 403406, 270.54861, 0.9287892); |
| | 0 | 282 | | sum += PeriodicTerm(julianCenturies, 195207, 340.19128, 35999.1376958); |
| | 0 | 283 | | sum += PeriodicTerm(julianCenturies, 119433, 63.91854, 35999.4089666); |
| | 0 | 284 | | sum += PeriodicTerm(julianCenturies, 112392, 331.2622, 35998.7287385); |
| | 0 | 285 | | sum += PeriodicTerm(julianCenturies, 3891, 317.843, 71998.20261); |
| | 0 | 286 | | sum += PeriodicTerm(julianCenturies, 2819, 86.631, 71998.4403); |
| | 0 | 287 | | sum += PeriodicTerm(julianCenturies, 1721, 240.052, 36000.35726); |
| | 0 | 288 | | sum += PeriodicTerm(julianCenturies, 660, 310.26, 71997.4812); |
| | 0 | 289 | | sum += PeriodicTerm(julianCenturies, 350, 247.23, 32964.4678); |
| | 0 | 290 | | sum += PeriodicTerm(julianCenturies, 334, 260.87, -19.441); |
| | 0 | 291 | | sum += PeriodicTerm(julianCenturies, 314, 297.82, 445267.1117); |
| | 0 | 292 | | sum += PeriodicTerm(julianCenturies, 268, 343.14, 45036.884); |
| | 0 | 293 | | sum += PeriodicTerm(julianCenturies, 242, 166.79, 3.1008); |
| | 0 | 294 | | sum += PeriodicTerm(julianCenturies, 234, 81.53, 22518.4434); |
| | 0 | 295 | | sum += PeriodicTerm(julianCenturies, 158, 3.5, -19.9739); |
| | 0 | 296 | | sum += PeriodicTerm(julianCenturies, 132, 132.75, 65928.9345); |
| | 0 | 297 | | sum += PeriodicTerm(julianCenturies, 129, 182.95, 9038.0293); |
| | 0 | 298 | | sum += PeriodicTerm(julianCenturies, 114, 162.03, 3034.7684); |
| | 0 | 299 | | sum += PeriodicTerm(julianCenturies, 99, 29.8, 33718.148); |
| | 0 | 300 | | sum += PeriodicTerm(julianCenturies, 93, 266.4, 3034.448); |
| | 0 | 301 | | sum += PeriodicTerm(julianCenturies, 86, 249.2, -2280.773); |
| | 0 | 302 | | sum += PeriodicTerm(julianCenturies, 78, 157.6, 29929.992); |
| | 0 | 303 | | sum += PeriodicTerm(julianCenturies, 72, 257.8, 31556.493); |
| | 0 | 304 | | sum += PeriodicTerm(julianCenturies, 68, 185.1, 149.588); |
| | 0 | 305 | | sum += PeriodicTerm(julianCenturies, 64, 69.9, 9037.75); |
| | 0 | 306 | | sum += PeriodicTerm(julianCenturies, 46, 8.0, 107997.405); |
| | 0 | 307 | | sum += PeriodicTerm(julianCenturies, 38, 197.1, -4444.176); |
| | 0 | 308 | | sum += PeriodicTerm(julianCenturies, 37, 250.4, 151.771); |
| | 0 | 309 | | sum += PeriodicTerm(julianCenturies, 32, 65.3, 67555.316); |
| | 0 | 310 | | sum += PeriodicTerm(julianCenturies, 29, 162.7, 31556.08); |
| | 0 | 311 | | sum += PeriodicTerm(julianCenturies, 28, 341.5, -4561.54); |
| | 0 | 312 | | sum += PeriodicTerm(julianCenturies, 27, 291.6, 107996.706); |
| | 0 | 313 | | sum += PeriodicTerm(julianCenturies, 27, 98.5, 1221.655); |
| | 0 | 314 | | sum += PeriodicTerm(julianCenturies, 25, 146.7, 62894.167); |
| | 0 | 315 | | sum += PeriodicTerm(julianCenturies, 24, 110.0, 31437.369); |
| | 0 | 316 | | sum += PeriodicTerm(julianCenturies, 21, 5.2, 14578.298); |
| | 0 | 317 | | sum += PeriodicTerm(julianCenturies, 21, 342.6, -31931.757); |
| | 0 | 318 | | sum += PeriodicTerm(julianCenturies, 20, 230.9, 34777.243); |
| | 0 | 319 | | sum += PeriodicTerm(julianCenturies, 18, 256.1, 1221.999); |
| | 0 | 320 | | sum += PeriodicTerm(julianCenturies, 17, 45.3, 62894.511); |
| | 0 | 321 | | sum += PeriodicTerm(julianCenturies, 14, 242.9, -4442.039); |
| | 0 | 322 | | sum += PeriodicTerm(julianCenturies, 13, 115.2, 107997.909); |
| | 0 | 323 | | sum += PeriodicTerm(julianCenturies, 13, 151.8, 119.066); |
| | 0 | 324 | | sum += PeriodicTerm(julianCenturies, 13, 285.3, 16859.071); |
| | 0 | 325 | | sum += PeriodicTerm(julianCenturies, 12, 53.3, -4.578); |
| | 0 | 326 | | sum += PeriodicTerm(julianCenturies, 10, 126.6, 26895.292); |
| | 0 | 327 | | sum += PeriodicTerm(julianCenturies, 10, 205.7, -39.127); |
| | 0 | 328 | | sum += PeriodicTerm(julianCenturies, 10, 85.9, 12297.536); |
| | 0 | 329 | | sum += PeriodicTerm(julianCenturies, 10, 146.1, 90073.778); |
| | 0 | 330 | | return sum; |
| | | 331 | | } |
| | | 332 | | |
| | | 333 | | private static double Aberration(double julianCenturies) |
| | | 334 | | { |
| | 0 | 335 | | return (0.0000974 * CosOfDegree(177.63 + (35999.01848 * julianCenturies))) - 0.005575; |
| | | 336 | | } |
| | | 337 | | |
| | | 338 | | private static double Nutation(double julianCenturies) |
| | | 339 | | { |
| | 0 | 340 | | double a = PolynomialSum(CoefficientsA, julianCenturies); |
| | 0 | 341 | | double b = PolynomialSum(CoefficientsB, julianCenturies); |
| | 0 | 342 | | return (-0.004778 * SinOfDegree(a)) - (0.0003667 * SinOfDegree(b)); |
| | | 343 | | } |
| | | 344 | | |
| | | 345 | | public static double Compute(double time) |
| | | 346 | | { |
| | 0 | 347 | | double julianCenturies = JulianCenturies(time); |
| | 0 | 348 | | double lambda = 282.7771834 |
| | 0 | 349 | | + (36000.76953744 * julianCenturies) |
| | 0 | 350 | | + (0.000005729577951308232 * SumLongSequenceOfPeriodicTerms(julianCenturies)); |
| | | 351 | | |
| | 0 | 352 | | double longitude = lambda + Aberration(julianCenturies) + Nutation(julianCenturies); |
| | 0 | 353 | | return InitLongitude(longitude); |
| | | 354 | | } |
| | | 355 | | |
| | | 356 | | public static double AsSeason(double longitude) |
| | | 357 | | { |
| | 0 | 358 | | return (longitude < 0) ? (longitude + FullCircleOfArc) : longitude; |
| | | 359 | | } |
| | | 360 | | |
| | | 361 | | private static double EstimatePrior(double longitude, double time) |
| | | 362 | | { |
| | 0 | 363 | | double timeSunLastAtLongitude = time - (MeanSpeedOfSun * AsSeason(InitLongitude(Compute(time) - longitude))) |
| | 0 | 364 | | double longitudeErrorDelta = InitLongitude(Compute(timeSunLastAtLongitude) - longitude); |
| | 0 | 365 | | return Math.Min(time, timeSunLastAtLongitude - (MeanSpeedOfSun * longitudeErrorDelta)); |
| | | 366 | | } |
| | | 367 | | |
| | | 368 | | // persian-new-year-on-or-before |
| | | 369 | | // number of days is the absolute date. The absolute date is the number of days from January 1st, 1 A.D. |
| | | 370 | | // 1/1/0001 is absolute date 1. |
| | | 371 | | internal static long PersianNewYearOnOrBefore(long numberOfDays) |
| | | 372 | | { |
| | 0 | 373 | | double date = (double)numberOfDays; |
| | | 374 | | |
| | 0 | 375 | | double approx = EstimatePrior(LongitudeSpring, MiddayAtPersianObservationSite(date)); |
| | 0 | 376 | | long lowerBoundNewYearDay = (long)Math.Floor(approx) - 1; |
| | 0 | 377 | | long upperBoundNewYearDay = lowerBoundNewYearDay + 3; // estimate is generally within a day of the actual oc |
| | 0 | 378 | | long day = lowerBoundNewYearDay; |
| | 0 | 379 | | for (; day != upperBoundNewYearDay; ++day) |
| | | 380 | | { |
| | 0 | 381 | | double midday = MiddayAtPersianObservationSite((double)day); |
| | 0 | 382 | | double l = Compute(midday); |
| | 0 | 383 | | if ((LongitudeSpring <= l) && (l <= TwoDegreesAfterSpring)) |
| | | 384 | | { |
| | | 385 | | break; |
| | | 386 | | } |
| | | 387 | | } |
| | 0 | 388 | | Debug.Assert(day != upperBoundNewYearDay); |
| | | 389 | | |
| | 0 | 390 | | return day - 1; |
| | | 391 | | } |
| | | 392 | | } |
| | | 393 | | } |
| | | 394 | | |