Geek Answers Handbook

Sunrise time in c

In my C application, I want to calculate sunrise / sunset times for a given date, latitude and longitude. I searched on the net, but I can not find a working sample.

I tried to implement this example: http://souptonuts.sourceforge.net/code/sunrise.c.html

But this sample did not work correctly.

Is there any simple source code or method C that can be easily implemented in my application?

Edit:
I am implementing the code on this link , but this gave me the wrong sunset / sunrise values. I also tried the Saul link here , but it also gave me the wrong result.
I have 41N, 28E. When I try to use the codes, both samples say that the sunrise is approximately 10:13 and the sunset is 23:24. But the correct values ​​are 06:06, 20:13.
I can not understand the problem.

+5
source share
11 answers

The sample code works in VC ++ 2010 with a few minor changes:

  • Compile it as a C ++ file, not C.
  • Delete the line #include <sys/time.h>.
  • Add #define _USE_MATH_DEFINESto the beginning of the file to determine the M_PI.
  • %T strftime() %X.

, , , , . printf(), .

, ( , , ).

+1

/ ,

  • N1 = (275 * /9) N2 = (( + 9)/12) N3 = (1 + (( - 4 * (/4) + 2)/3)) N = N1 - (N2 * N3) + - 30

  • lngHour = /15

    : t = N + ((6 - lngHour)/24) : t = N + ((18 - lngHour)/24)

  • M = (0,9856 * t) - 3,259

  • L = M + (1,916 * sin (M)) + (0,020 * sin (2 * M)) + 282,634 : L [0,360), / 360

5.

RA = atan(0.91764 * tan(L))
NOTE: RA potentially needs to be adjusted into the range [0,360) by adding/subtracting 360

5. , L

Lquadrant  = (floor( L/90)) * 90
RAquadrant = (floor(RA/90)) * 90
RA = RA + (Lquadrant - RAquadrant)

5.

RA = RA / 15

  1. sinDec = 0,399782 * sin (L) cosDec = cos (asin (sinDec))

7.

cosH = (cos(zenith) - (sinDec * sin(latitude))) / (cosDec * cos(latitude))

if (cosH >  1) 
  the sun never rises on this location (on the specified date)
if (cosH < -1)
  the sun never sets on this location (on the specified date)

7b. H

if if rising time is desired:
  H = 360 - acos(cosH)
if setting time is desired:
  H = acos(cosH)

H = H / 15

  1. /

    T = H + RA - (0,06571 * t) - 6,622

  2. UTC

    UT = T - lngHour : UT [0,24] / 24

  3. UT /

    localt= UT + localOffset

+8

, . . , ...

#include "stdafx.h"
#include <iostream>
#include <math.h>
#include <time.h>

using namespace std;


//STANDARD CONSTANTS
double pi = 3.1415926535;   // Pi
double solarConst = 1367;           // solar constant W.m-2


// Function to convert radian to hours
double RadToHours (double tmp)
{
    //double pi = 3.1415926535; // Pi
    return (tmp * 12 / pi);
}
// Function to convert hours to radians
double HoursToRads (double tmp)
{
    //double pi = 3.1415926535; // Pi
    return (tmp * pi / 12);
}




// Function to calculate the angle of the day
double AngleOfDay (int day,     // number of the day 
                   int month,   // number of the month
                   int year // year 
                 )

{   // local vars
    int i, leap;
    int numOfDays = 0;                                              // number of Day 13 Nov=317
    int numOfDaysofMonths[12] = {0,31,28,31,30,31,30,31,31,30,31,30};   // Number of days per month
    int AllYearDays;                                                // Total number of days in a year 365 or 366
    double DayAngle;                                                // angle of the day (radian)
    //double pi = 3.1415926535; // Pi

    // leap year ??
    leap = 0;
    if ((year % 400)==0) 
    {   AllYearDays = 366;
        leap = 1;
    }
    else if ((year % 100)==0) AllYearDays = 365;
         else if ((year % 4)==0)
            {   AllYearDays = 366;
                leap = 1;
            }
             else AllYearDays = 365;

    // calculate number of day
    for (i=0;i<month;i++) numOfDays += numOfDaysofMonths[i];
    if ( (month > 2) && leap) numOfDays++;
    numOfDays += day;

    // calculate angle of day
    DayAngle = (2*pi*(numOfDays-1)) / AllYearDays;
    return DayAngle;

}


// Function to calculate declination - in radian
double Declination (double DayAngle     // angle day in radian
                    )
{
    double SolarDeclination;
    // Solar declination (radian)
    SolarDeclination = 0.006918 
            - 0.399912 * cos (DayAngle)
            + 0.070257 * sin (DayAngle)
            - 0.006758 * cos (2*DayAngle)
            + 0.000907 * sin (2*DayAngle)
            - 0.002697 * cos (3*DayAngle)
            + 0.00148 * sin (3*DayAngle);
    return SolarDeclination;
}



// Function to calculate Equation of time ( et = TSV - TU )
double EqOfTime (double DayAngle        // angle day (radian)
                      )
{
    double et;
    // Equation of time (radian)
    et = 0.000075
         + 0.001868 * cos (DayAngle)
         - 0.032077 * sin (DayAngle)
         - 0.014615 * cos (2*DayAngle)
         - 0.04089 * sin (2*DayAngle);
    // Equation of time in hours
    et = RadToHours(et);

    return et;
}

// Calculation of the duration of the day in radian
double DayDurationRadian (double _declination,      // _declination in radian
                          double lat                // latitude in radian
                         )
{
    double dayDurationj;

    dayDurationj = 2 * acos( -tan(lat) * tan(_declination) );
    return dayDurationj;
}

// Function to calculate Day duration in Hours
double DayDuratInHours (double _declination     // _declination in radian
                  , double lat              // latitude in radian
                  )
{
    double dayDurationj;

    dayDurationj = DayDurationRadian(_declination, lat);
    dayDurationj = RadToHours(dayDurationj);
    return dayDurationj;
}


// Function to calculate the times TSV-UTC
double Tsv_Tu (double rlong             // longitude en radian positive a l est. 
               ,double eqOfTime         // Equation of times en heure
              )
{
    double diffUTC_TSV; double pi = 3.1415926535;   // Pi

    // diffUTC_TSV Solar time as a function of longitude and the eqation of time
    diffUTC_TSV = rlong * (12 / pi) + eqOfTime;

    // difference with local time
    return diffUTC_TSV;
}


// Calculations of the orbital excentricity
double Excentricity(int day,
                    int month,
                    int year)
{

    double dayAngleRad, E0;

    // calculate the angle of day in radian
    dayAngleRad = AngleOfDay(day, month, year);

    // calculate the excentricity
    E0 = 1.000110 + 0.034221 * cos(dayAngleRad) 
            + 0.001280 * sin(dayAngleRad)
            +0.000719 * cos(2*dayAngleRad)
            +0.000077 * sin(2*dayAngleRad);

    return E0;
}

// Calculate the theoretical energy flux for the day radiation
double TheoreticRadiation(int day, int month, int year, 
                            double lat          // Latitude in radian !
                            )
{
    double RGth;        // Theoretical radiation
    double decli;       // Declination
    double E0;
    double sunriseHourAngle;            // Hour angle of sunset



    // Calculation of the declination in radian
    decli = Declination (AngleOfDay(day, month, year));

    // Calcuate excentricity
    E0 = Excentricity(day, month, year);

    // Calculate hour angle in radian
    sunriseHourAngle = DayDurationRadian(decli, lat) / 2;

    // Calculate Theoretical radiation en W.m-2
    RGth = solarConst * E0 * (cos(decli)*cos(lat)*sin(sunriseHourAngle)/sunriseHourAngle + sin(decli)*sin(lat));

    return RGth;

}

// Function to calculate decimal hour of sunrise: result in local hour
double CalclulateSunriseLocalTime(int day,
                        int month,
                        int year,
                        double rlong,
                        double rlat)

{   
    // local variables
    int h1, h2;
    time_t hour_machine;
    struct tm *local_hour, *gmt_hour;


    double result;
    // Calculate the angle of the day
    double DayAngle = AngleOfDay(day, month, year);
    // Declination
    double SolarDeclination = Declination(DayAngle);
    // Equation of times
    double eth = EqOfTime(DayAngle);
    // True solar time
    double diffUTC_TSV = Tsv_Tu(rlong,eth);
    // Day duration
    double dayDurationj = DayDuratInHours(SolarDeclination,rlat);

    // local time adjust
    time( &hour_machine );      // Get time as long integer. 
    gmt_hour =  gmtime( &hour_machine );
    h1 = gmt_hour->tm_hour;
    local_hour = localtime( &hour_machine );    // local time. 
    h2 = local_hour->tm_hour;

    // final result
    result = 12 - fabs(dayDurationj / 2) - diffUTC_TSV + h2-h1;

    return result;

}

// Function to calculate decimal hour of sunset: result in local hour
double CalculateSunsetLocalTime(int day,
                          int month,
                          int year,
                          double rlong,
                          double rlat)

{   
    // local variables
    int h1, h2;
    time_t hour_machine;
    struct tm *local_hour, *gmt_hour;
    double result;

    // Calculate the angle of the day
    double DayAngle = AngleOfDay(day, month, year);
    // Declination
    double SolarDeclination = Declination(DayAngle);
    // Equation of times
    double eth = EqOfTime(DayAngle);
    // True solar time
    double diffUTC_TSV = Tsv_Tu(rlong,eth);
    // Day duration
    double dayDurationj = DayDuratInHours(SolarDeclination,rlat);

    // local time adjust
    time( &hour_machine );      // Get time as long integer. 
    gmt_hour =  gmtime( &hour_machine );
    h1 = gmt_hour->tm_hour;
    local_hour = localtime( &hour_machine );    // local time. 
    h2 = local_hour->tm_hour;

    // resultat
    result = 12 + fabs(dayDurationj / 2) - diffUTC_TSV + h2-h1;

    return result;

}



// Function to calculate decimal hour of sunrise: result universal time
double CalculateSunriseUniversalTime(int day,
                        int month,
                        int year,
                        double rlong,
                        double rlat)

{   
    double result;
    // Calculate the angle of the day
    double DayAngle = AngleOfDay(day, month, year);
    // Declination
    double SolarDeclination = Declination(DayAngle);
    // Equation of times
    double eth = EqOfTime(DayAngle);
    // True solar time
    double diffUTC_TSV = Tsv_Tu(rlong,eth);
    // Day duration
    double dayDurationj = DayDuratInHours(SolarDeclination,rlat);
    // resultat
    result = 12 - fabs(dayDurationj / 2) - diffUTC_TSV;

    return result;

}

// Function to calculate decimal hour of sunset: result in universal time
double CalculateSunsetUniversalTime(int day,
                          int month,
                          int year,
                          double rlong,
                          double rlat)

{   
    double result;

    // Calculate the angle of the day
    double DayAngle = AngleOfDay(day, month, year);
    // Declination
    double SolarDeclination = Declination(DayAngle);
    // Equation of times
    double eth = EqOfTime(DayAngle);
    // True solar time
    double diffUTC_TSV = Tsv_Tu(rlong,eth);
    // Day duration
    double dayDurationj = DayDuratInHours(SolarDeclination,rlat);
    // resultat
    result = 12 + fabs(dayDurationj / 2) - diffUTC_TSV;

    return result;

}


// Function to calculate the height of the sun in radians the day to day j and hour TU
double SolarHeight (int tu,     // universal times (0,1,2,.....,23)
                      int day,
                      int month,
                      int year,
                      double lat,   // latitude in radian
                      double rlong  // longitude in radian
                     )
{
    // local variables
    double pi = 3.1415926535;   // Pi
    double result, tsvh;

    // angle of the day
    double DayAngle = AngleOfDay(day, month, year);
    // _declination
    double decli = Declination(DayAngle);
    // eq of time
    double eq = EqOfTime(DayAngle);
    // calculate the tsvh with rlong positiv for the east and negative for the west
    tsvh = tu + rlong*180/(15*pi) + eq;
    // hour angle per hour
    double ah = acos( -cos((pi/12)*tsvh) );
    // final result
    result = asin( sin(lat)*sin(decli) + cos(lat)*cos(decli)*cos(ah) );


    return result;
}



///////////EXTRA FUNCTIONS/////////////////////////////

//Julian day conversion for days calculations
//Explanation for this sick formula...for the curious guys...
//http://www.cs.utsa.edu/~cs1063/projects/Spring2011/Project1/jdn-explanation.html

int julian(int year, int month, int day) {
  int a = (14 - month) / 12;
  int y = year + 4800 - a;
  int m = month + 12 * a - 3;
  if (year > 1582 || (year == 1582 && month > 10) || (year == 1582 && month == 10 && day >= 15))
    return day + (153 * m + 2) / 5 + 365 * y + y / 4 - y / 100 + y / 400 - 32045;
  else
    return day + (153 * m + 2) / 5 + 365 * y + y / 4 - 32083;
}



int _tmain(int argc, _TCHAR* argv[])
{
    int day = 14;
    int month = 11;
    int year = 2013;
    double lat = 39.38;
    double lon = 22.75;
    double rlat = 39.38 * pi/180;
    double rlong = 22.75 * pi/180;

    double _AngleOfDay = AngleOfDay ( day , month , year );
    cout << "Angle of day: " << _AngleOfDay << "\n";

    double _Declinaison = Declination (_AngleOfDay);
    cout << "Declination (Delta): " << _Declinaison << "\n";

    double _EqOfTime = EqOfTime (_AngleOfDay);
    cout << "Declination (Delta): " << _EqOfTime << "\n";

    double _DayDuratInHours = DayDuratInHours (_Declinaison, rlat);
    cout << "Day duration: " << _DayDuratInHours << "\n";

    double _Excentricity = Excentricity(day, month, year);
    cout << "Excentricity: " << _Excentricity << "\n";

    double _TheoreticRadiation = TheoreticRadiation(day, month, year, rlat);
    cout << "Theoretical radiation: " << _TheoreticRadiation << "\n";

    double _CalclulateSunriseLocalTime = CalclulateSunriseLocalTime
                    (day, month, year, rlong, rlat);
    cout << "Sunrise Local Time: " << _CalclulateSunriseLocalTime << "\n"; 

    double _CalculateSunsetLocalTime = CalculateSunsetLocalTime
        (day, month, year, rlong, rlat);
    cout << "Sunrise Local Time: " << _CalculateSunsetLocalTime << "\n";

    return 0;
}
+2

C, -, , C++, C#, :

+1

guide ( @BenjaminMonate, , @Geetha ), C , , .

#include <math.h>
#define PI 3.1415926
#define ZENITH -.83

float calculateSunrise(int year,int month,int day,float lat, float lng,int localOffset, int daylightSavings) {
    /*
    localOffset will be <0 for western hemisphere and >0 for eastern hemisphere
    daylightSavings should be 1 if it is in effect during the summer otherwise it should be 0
    */
    //1. first calculate the day of the year
    float N1 = floor(275 * month / 9);
    float N2 = floor((month + 9) / 12);
    float N3 = (1 + floor((year - 4 * floor(year / 4) + 2) / 3));
    float N = N1 - (N2 * N3) + day - 30;

    //2. convert the longitude to hour value and calculate an approximate time
    float lngHour = lng / 15.0;      
    float t = N + ((6 - lngHour) / 24);   //if rising time is desired:
    //float t = N + ((18 - lngHour) / 24)   //if setting time is desired:

    //3. calculate the Sun mean anomaly   
    float M = (0.9856 * t) - 3.289;

    //4. calculate the Sun true longitude
    float L = fmod(M + (1.916 * sin((PI/180)*M)) + (0.020 * sin(2 *(PI/180) * M)) + 282.634,360.0);

    //5a. calculate the Sun right ascension      
    float RA = fmod(180/PI*atan(0.91764 * tan((PI/180)*L)),360.0);

    //5b. right ascension value needs to be in the same quadrant as L   
    float Lquadrant  = floor( L/90) * 90;
    float RAquadrant = floor(RA/90) * 90;
    RA = RA + (Lquadrant - RAquadrant);

    //5c. right ascension value needs to be converted into hours   
    RA = RA / 15;

    //6. calculate the Sun declination
    float sinDec = 0.39782 * sin((PI/180)*L);
    float cosDec = cos(asin(sinDec));

    //7a. calculate the Sun local hour angle
    float cosH = (sin((PI/180)*ZENITH) - (sinDec * sin((PI/180)*lat))) / (cosDec * cos((PI/180)*lat));
    /*   
    if (cosH >  1) 
    the sun never rises on this location (on the specified date)
    if (cosH < -1)
    the sun never sets on this location (on the specified date)
    */

    //7b. finish calculating H and convert into hours
    float H = 360 - (180/PI)*acos(cosH);   //   if if rising time is desired:
    //float H = acos(cosH) //   if setting time is desired:      
    H = H / 15;

    //8. calculate local mean time of rising/setting      
    float T = H + RA - (0.06571 * t) - 6.622;

    //9. adjust back to UTC
    float UT = fmod(T - lngHour,24.0);

    //10. convert UT value to local time zone of latitude/longitude
    return UT + localOffset + daylightSavings;

    }

void printSunrise() {
    float localT = calculateSunrise(/*args*/);
    double hours;
    float minutes = modf(localT,&hours)*60;
    printf("%.0f:%.0f",hours,minutes);
    }
+1

, scottmrogowski, ,

  • //float H = (180/PI) * acos (cosH)// :
  • float localt= fmod (24 + calculateSunrise (/* args */), 24.0);// printSunrise

Marek

+1

Tomoyose C #. , - 4 -. , , , - . , - - , . - , : -)

NodaTime, Java Skeet Java JodaTime. NodaTime nuGet.

using System;

using NodaTime;

namespace YourNamespaceHere
{
    public static class MySunset
    {
        // Based on Tomoyose's
        // http://stackoverflow.com/questions/7064531/sunrise-sunset-times-in-c

        private static DateTimeZone mUtcZone = DateTimeZoneProviders.Tzdb["Etc/UTC"];
        private static int mSecondsInDay = 24 * 60 * 60;

        // Convert radian to hours
        private static double RadiansToHours(double radians)
        {
            return (radians * 12.0 / Math.PI);
        }

        // Convert hours to radians
        private static double HoursToRadians(double hours)
        {
            return (hours * Math.PI / 12.0);
        }

        // Calculate the angle of the day
        private static double AngleOfDay(int year, int month, int day) 
        {   
            DateTime date = new DateTime(year, month, day);
            int daysInYear = DateTime.IsLeapYear(year) ? 366 : 365;

            // Return angle of day in radians
            return (2.0 * Math.PI * ((double)date.DayOfYear - 1.0)) / (double)daysInYear;
        }

        // Calculate declination in radians
        private static double SolarDeclination(double angleOfDayRadians)
        {
            // Return solar declination in radians
            return 0.006918
                - 0.399912 * Math.Cos(angleOfDayRadians)
                + 0.070257 * Math.Sin(angleOfDayRadians)
                - 0.006758 * Math.Cos(2.0 * angleOfDayRadians)
                + 0.000907 * Math.Sin(2.0 * angleOfDayRadians)
                - 0.002697 * Math.Cos(3.0 * angleOfDayRadians)
                + 0.00148 * Math.Sin(3.0 * angleOfDayRadians);
        }

        // Calculate Equation of time ( eot = TSV - TU )
        private static double EquationOfTime(double angleOfDayRadians)
        {
            // Equation of time (radians)
            double et = 0.000075
                 + 0.001868 * Math.Cos(angleOfDayRadians)
                 - 0.032077 * Math.Sin(angleOfDayRadians)
                 - 0.014615 * Math.Cos(2.0 * angleOfDayRadians)
                 - 0.04089 * Math.Sin(2.0 * angleOfDayRadians);

            // Return equation-of-time in hours
            return RadiansToHours(et);
        }

        // Calculate the duration of the day in radians
        private static double DayDurationRadians(double declinationRadians, double latitudeRadians)
        {
            return 2.0 * Math.Acos(-Math.Tan(latitudeRadians) * Math.Tan(declinationRadians));
        }

        // Calculate day duration in hours
        private static double DayDurationHours(double declinationRadians, double latitudeRadians)
        {
            return RadiansToHours(
                DayDurationRadians(declinationRadians, latitudeRadians)
            );
        }

        // Calculate the times TSV-UTC
        private static double Tsv_Tu(double longitudeRadians, double equationOfTime)
        {
            // Solar time as a function of longitude and the equation of time
            return longitudeRadians * (12.0 / Math.PI) + equationOfTime;
        }

        private static void GetDayParameters(int year, int month, int day, double latitude, double longitude,
            out double dayDuration, out double diffUTC_TSV)
        {
            double latitudeRadians = latitude * Math.PI / 180.0;
            double longitudeRadians = longitude * Math.PI / 180.0;

            // Calculate the angle of the day
            double dayAngle = AngleOfDay(year, month, day);
            // Declination
            double solarDeclination = SolarDeclination(dayAngle);
            // Equation of times
            double equationOfTime = EquationOfTime(dayAngle);
            // True solar time
            diffUTC_TSV = Tsv_Tu(longitudeRadians, equationOfTime);
            // Day duration
            dayDuration = DayDurationHours(solarDeclination, latitudeRadians);
        }

        // Calculate decimal UTC hour of sunrise.
        private static double CalculateSunriseUTC(int year, int month, int day, double latitude, double longitude)
        {
            double dayDuration;
            double diffUTC_TSV;

            GetDayParameters(year, month, day, latitude, longitude, out dayDuration, out diffUTC_TSV);

            return 12.0 - Math.Abs(dayDuration / 2.0) - diffUTC_TSV;
        }

        // Calculate decimal UTC hour of sunset. 
        private static double CalculateSunsetUTC(int year, int month, int day, double latitude, double longitude)
        {
            double dayDuration;
            double diffUTC_TSV;

            GetDayParameters(year, month, day, latitude, longitude, out dayDuration, out diffUTC_TSV);

            return 12.0 + Math.Abs(dayDuration / 2.0) - diffUTC_TSV;
        }

        // Public methods to return sun rise and set times.
        public static Tuple<ZonedDateTime, ZonedDateTime> GetSunRiseSet(LocalDate dateAtLocation, DateTimeZone locationZone, double latitude, double longitude)
        {
            // latitude-longitude must lie within zone of locationZone


            // Get UTC rise and set
            double dayDuration;
            double diffUTC_TSV;
            GetDayParameters(dateAtLocation.Year, dateAtLocation.Month, dateAtLocation.Day, latitude, longitude, out dayDuration, out diffUTC_TSV);
            double sunriseUtcDecimal = 12.0 - Math.Abs(dayDuration / 2.0) - diffUTC_TSV;
            double sunsetUtcDecimal = 12.0 + Math.Abs(dayDuration / 2.0) - diffUTC_TSV;

            // Convert decimal UTC to UTC dates

            // If a UTC time is negative then it means the date before in the UTC timezone.
            // So if negative need to minus 1 day from date and set time as 24 - decimal_time.
            LocalDateTime utcRiseLocal;
            LocalDateTime utcSetLocal;
            if (sunriseUtcDecimal < 0)
            {
                LocalDate utcDateAdjusted = dateAtLocation.PlusDays(-1);
                // Normalize() is important here; otherwise only have access to seconds.
                Period utcTimeAdjusted = Period.FromSeconds((long)((24.0 + sunriseUtcDecimal) / 24.0 * mSecondsInDay)).Normalize(); // + a negative
                utcRiseLocal = new LocalDateTime(utcDateAdjusted.Year, utcDateAdjusted.Month, utcDateAdjusted.Day,
                    (int)utcTimeAdjusted.Hours, (int)utcTimeAdjusted.Minutes, (int)utcTimeAdjusted.Seconds);
            }
            else
            {
                Period utcTime = Period.FromSeconds((long)(sunriseUtcDecimal / 24.0 * mSecondsInDay)).Normalize();
                utcRiseLocal = new LocalDateTime(dateAtLocation.Year, dateAtLocation.Month, dateAtLocation.Day,
                    (int)utcTime.Hours, (int)utcTime.Minutes, (int)utcTime.Seconds);
            }
            if (sunsetUtcDecimal < 0) // Maybe not possible?
            {
                LocalDate utcDateAdjusted = dateAtLocation.PlusDays(-1);
                Period utcTimeAdjusted = Period.FromSeconds((long)((24.0 + sunsetUtcDecimal) / 24.0 * mSecondsInDay)).Normalize();
                utcSetLocal = new LocalDateTime(utcDateAdjusted.Year, utcDateAdjusted.Month, utcDateAdjusted.Day,
                    (int)utcTimeAdjusted.Hours, (int)utcTimeAdjusted.Minutes, (int)utcTimeAdjusted.Seconds);
            }
            else
            {
                Period utcTime = Period.FromSeconds((long)(sunsetUtcDecimal / 24.0 * mSecondsInDay)).Normalize();
                utcSetLocal = new LocalDateTime(dateAtLocation.Year, dateAtLocation.Month, dateAtLocation.Day,
                    (int)utcTime.Hours, (int)utcTime.Minutes, (int)utcTime.Seconds);
            }

            // FIX: always about 4 minutes later/earlier than other sources
            utcRiseLocal = utcRiseLocal.PlusMinutes(-4);
            utcSetLocal = utcSetLocal.PlusMinutes(4);

            // Get zoned datetime from UTC local datetimes
            ZonedDateTime utcRiseZoned = new LocalDateTime(utcRiseLocal.Year, utcRiseLocal.Month, utcRiseLocal.Day, utcRiseLocal.Hour, utcRiseLocal.Minute, utcRiseLocal.Second).InZoneLeniently(mUtcZone);
            ZonedDateTime utcSetZoned = new LocalDateTime(utcSetLocal.Year, utcSetLocal.Month, utcSetLocal.Day, utcSetLocal.Hour, utcSetLocal.Minute, utcSetLocal.Second).InZoneLeniently(mUtcZone);

            // Return zoned UTC to zoned local 
            return new Tuple<ZonedDateTime, ZonedDateTime>
            (
                new ZonedDateTime(utcRiseZoned.ToInstant(), locationZone),
                new ZonedDateTime(utcSetZoned.ToInstant(), locationZone)
            );
        }
        public static Tuple<ZonedDateTime, ZonedDateTime> GetSunRiseSet(int year, int month, int day, string locationZone, double latitude, double longitude)
        {
            return GetSunRiseSet(new LocalDate(year, month, day), DateTimeZoneProviders.Tzdb[locationZone], latitude, longitude);
        }
        public static Tuple<ZonedDateTime, ZonedDateTime> GetSunRiseSet(ZonedDateTime zonedDateAtLocation, double latitude, double longitude)
        {
            return GetSunRiseSet(zonedDateAtLocation.LocalDateTime.Date, zonedDateAtLocation.Zone, latitude, longitude);
        }
    }
}
+1

It is adaptable and pretty accurate. You will get all the components, and then all you need to calculate is the cosine of the arc for the Zenith to Horizon angle. Yes, there are simpler ways, but don't you basically want to track the sun? It may come in handy someday. http://www.nrel.gov/midc/spa/

0
source 
//a practical approximation for a microcontroller using lookup table
//seems to only use a small fraction of the resources required by the trigonometric functions
//ignores daylight saving: idea is for the device to approximate and trigger actual sunrise, sunset event as opposed to actual (politically correct local hhmm)
//using globals gps.Lat gps.LatDir(0:S 1:N) gps.Day (day of month) gps.Mth
//needs solar noon offset,latitude,dayOfMonth,Month
//LocalNoon-SolarNoon (may be off +/- up to ?? 2 hours) depending on local timezone and longitude (? typical discrepancy about 20min)

void SunriseSunsetCalculations(u8 SolarNoonOffsetInDeciHours,u8 SolarNoonOffsetDirection){
    if(SolarNoonOffsetInDeciHours>20){SolarNoonOffsetInDeciHours=0;}//limit offest to 2 hours: discard offset on error
    //--------references:
    //https://www.orchidculture.com/COD/daylength.html
    //http://aa.usno.navy.mil/data/docs/Dur_OneYear.php
    //SolarTime is up two two hours off in either direction depending on timezone:
    //http://blog.poormansmath.net/how-much-is-time-wrong-around-the-world/
    //------------------
    //lookUpTable Of daylength in decihours(6min)(fits well in u8) - 15day And 5 DegLat resolution
    //SolarNoonOffsetDirection: 0: negative(b4 local noon), 1: +ve (solar noon after local noon)
    u8 a[13][12]={//length of day in decihours
            //   Dec      Nov     Oct     Sep     Aug     Jul
            //Jan     Feb     Mar     Apr     May     June
            {120,120,120,120,120,120,120,120,120,120,120,120},  //lat  0____0
            {126,125,124,123,122,121,119,118,116,115,115,114},  //lat 10____1
            {129,128,127,125,123,121,119,117,115,113,112,111},  //lat 15____2
            {132,131,129,127,124,121,118,115,113,110,109,108},  //lat 20____3
            {135,134,131,128,125,122,118,114,111,108,106,105},  //lat 25____4
            {139,137,134,130,127,122,117,113,108,105,102,101},  //lat 30____5
            {143,141,137,133,128,123,117,111,106,102, 98, 97},  //lat 35____6
            {148,145,141,135,130,123,116,109,103, 98, 94, 92},  //lat 40____7
            {154,150,145,138,132,124,115,107,100, 93, 88, 86},  //lat 45____8
            {161,157,150,142,134,124,114,105, 96, 88, 82, 79},  //lat 50____9
            {170,165,156,146,137,125,113,102, 91, 81, 73, 69},  //lat 55___10
            {183,176,165,152,140,126,112, 98, 84, 72, 61, 56},  //lat 60___11
            {200,185,171,152,134,121,101, 84, 65, 53, 40, 33}   //lat 65___12
    };
    u8 b[]={6,12,17,22,27,32,37,42,47,52,57,62,90}; // latitude limit cutoffs to get index of lookUpTable
    u8 lat=gps.Lat/10000000; //lat stored in u32 to 7 decimals resolution (32bit unsigned integer)
    u8 i=0; while(b[i]<lat){i++;}  //get row index for daylength table
    u8 k,ix; //k: 15 day offset;    ix: column index for daylength table
    k=gps.Day/15;if(k){k=1;}//which half of the month (avoid k=2 eg:31/15)
    if(!gps.LatDir){ //0:southern latitudes
        if(gps.Mth<7){
            ix=(gps.Mth-1)*2+k;     //2 fields per month (k to select)
        }else{                      //beyond june, use row in reverse
            ix=11-(gps.Mth-7)*2-k;  //2 fields per month (k to select)
        }
    }else{           //1:northern latitudes
        if(gps.Mth<7){              //with first six month read rows in reverse
            ix=11-(gps.Mth-1)*2-k;  //2 fields per month (k to select)
        }else{                      //beyond june, use same row forward
            ix=(gps.Mth-7)*2+k;     //2 fields per month (k to select)
        }
    }
    //debug only:...dcI("\r\ni", i ,Blue,Red,1);dcI("ix", ix ,Blue,Red,1);dcI("a[i][ix]", a[i][ix] ,Blue,Red,1);
    u8 h=a[i][ix]/2; //HalfTheDayLightLength in deciHours
    u8 j[]={-1,1};   //multiplier: 0:(-)  1:(+)
    u8 sn=120+SolarNoonOffsetInDeciHours*j[SolarNoonOffsetDirection]; //Solar noon
    u8 srdh=sn-h;    //sunrise in deciHours = solarNoon minus HalfTheDayLightLength
    u8 ssdh=sn+h;    //sunset  in deciHours = solarNoon  plus HalfTheDayLightLength
    //debug only:...dcI("\r\nSunRiseDeciHours", srdh ,Blue,Red,1);dcI("SunSetDeciHours", ssdh ,Blue,Red,1);
    gps.HmSunRise=deciHourTohhmm(srdh);
    gps.HmSunSet =deciHourTohhmm(ssdh);
}
u16 deciHourTohhmm(u8 dh){ //return unsigned integer from 0 to 2400
   u16 h=(dh/10)*100; //hours: hh00
   u8  r= dh%10;      //fraction hour remainder to be converted to minutes
   u8  m= 6*r;        //60*r/10
   return(h+m);
}


/*
Example Output: (!!! solarNoonOffset kept at 0(ignored) for the below example)

:_(08474300)___:_(A)___:_(381234567)___:_(S)___:_(1431234567)___:_(E)___
:_(GPS OK)___:_(Sat)___:_(12)___:_(Aug)___:_(2017)___hhmm:_(847)___
//...........
i:_(7)___ix:_(9)___a[i][ix]:_(98)___
SunRiseDeciHours:_(71)___SunSetDeciHours:_(169)___
HmSunRise:_(706)___
HmSunSet:_(1654)___
//............same as above but change LatDir to 1 (northern hemisphere)
i:_(7)___ix:_(2)___a[i][ix]:_(141)___
SunRiseDeciHours:_(50)___SunSetDeciHours:_(190)___
HmSunRise:_(500)___
HmSunSet:_(1900)___
..........ignore dcI(...) it just a custom function printing through the serial port
*/
0
source

More articles:


All Articles