SunriseSunset.m

BOOL isCurrentlyDaytime(CGFloat latitude, CGFloat longitude)
{
    // http://users.electromagnetic.net/bu/astro/iyf-calc.php
    // http://aa.quae.nl/en/reken/zonpositie.html
    
    static NSUInteger const secondsPerDay = 60*60*24;
    static NSTimeInterval const julianLeap = 0.0009;
    static double const degreeToRadian = M_PI/180.0;
    static double const radianToDegree = 180.0/M_PI;
    
    //reference date is Jan 1, 2000
    static NSUInteger const julianReferenceDate = 2451545;
    static NSTimeInterval const julianReferenceIntervalSince1970 = 946684800;
    
    static double const earthMeanAnomaly = 357.5291;
    static double const earthMeanAnomalyPerDay = 0.98560028;
    static double const earthEquationOfCenterCoefficient1 = 1.9148;
    static double const earthEquationOfCenterCoefficient2 = 0.0200;
    static double const earthEquationOfCenterCoefficient3 = 0.0003;
    static double const earthEccentricityVariance = 0.0053;
    static double const earthObliquityVariance = -0.0069;
    static double const earthObliquityEquator = 23.45;
    static double const earthEllipticLongitude = 102.9372;
    static double const earthLightRefraction = -0.83;
    
    NSDate *referenceDate = [NSDate dateWithTimeIntervalSince1970:julianReferenceIntervalSince1970];
    
    //julian cycle since Jan 1, 2000
    NSTimeInterval intervalSinceReferenceDate = [[NSDate date] timeIntervalSinceDate:referenceDate];
    NSUInteger daysSinceReferenceDate = intervalSinceReferenceDate / secondsPerDay;
    NSUInteger n = (NSUInteger)round((daysSinceReferenceDate-julianLeap) - (longitude / 360.0));
    
    double j1 = julianReferenceDate+julianLeap+(longitude/360.0)+n;
    double M = (earthMeanAnomaly + earthMeanAnomalyPerDay * (j1 - julianReferenceDate));
    M = fmod(M, 360.0);
    
    double C = (earthEquationOfCenterCoefficient1 * sin(degreeToRadian * M))
                        + (earthEquationOfCenterCoefficient2 * sin(degreeToRadian * 2.0 * M))
                        + (earthEquationOfCenterCoefficient3 * sin(degreeToRadian * 3.0 * M));
    double ellipticalLongitude = (M + earthEllipticLongitude + C + 180.0);
    ellipticalLongitude = fmod(ellipticalLongitude, 360.0);
    
    double julianTransit = j1
                            + (earthEccentricityVariance * sin(degreeToRadian * M))
                            + (earthObliquityVariance * sin(degreeToRadian * 2.0 * ellipticalLongitude));
    double solarDeclination = radianToDegree*asin(sin(degreeToRadian * ellipticalLongitude)
                                                      * sin(degreeToRadian * earthObliquityEquator));
    double hourAngle = radianToDegree*acos((sin(degreeToRadian * earthLightRefraction) - sin(degreeToRadian * latitude) * sin(degreeToRadian * solarDeclination))
                                               / (cos(degreeToRadian * latitude) * cos(degreeToRadian * solarDeclination)));
    double j2 = julianReferenceDate + julianLeap + ((hourAngle + longitude)/360.0) + n;
    
    NSTimeInterval julianSunset = j2
                                    + (earthEccentricityVariance * sin(degreeToRadian * M))
                                    + (earthObliquityVariance * sin(degreeToRadian * 2.0 * ellipticalLongitude));
    NSTimeInterval julianSunrise = julianTransit - (julianSunset - julianTransit);
    
    // http://stackoverflow.com/a/27709317/3813982
    NSTimeInterval JD_JAN_1_2000_0000GMT = julianReferenceDate - 0.5;
    NSTimeInterval sunriseInterval = (julianSunrise - JD_JAN_1_2000_0000GMT) * secondsPerDay;
    NSTimeInterval sunsetInterval = (julianSunset - JD_JAN_1_2000_0000GMT) * secondsPerDay;
    
#if DEBUG
    NSDate *sunrise = [NSDate dateWithTimeIntervalSince1970:julianReferenceIntervalSince1970+sunriseInterval];
    NSDate *sunset = [NSDate dateWithTimeIntervalSince1970:julianReferenceIntervalSince1970+sunsetInterval];
    NSDateFormatter *dateFormatter = [[NSDateFormatter alloc] init];
    [dateFormatter setDateFormat:@"hh:mm a"];
    __unused NSString *formattedSunrise = [dateFormatter stringFromDate:sunrise];
    __unused NSString *formattedSunset = [dateFormatter stringFromDate:sunset];
#endif
    
    return sunriseInterval <= intervalSinceReferenceDate && intervalSinceReferenceDate < sunsetInterval;
}