19315321183

Committed 12 Nov 2025 11:35PM UTC coverage: 88.752% (-1.3%) from 90.052%

Build # 19315321183

Build Type

Pull #83

github

Committed by

web-flow

Commit Message

Merge bfed87b33 into e67f9eadc

Pull Request Pull Request #83: Sun position calcs

Coverage Stats

996 of 1195 new or added lines in 4 files covered. (83.35%)

5492 of 6188 relevant lines covered (88.75%)

8572678.39 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

94.44

/coretrace/simulation_data/solar_position_calculators/basic_sun_position.cpp


#include "basic_sun_position.hpp"
#include "constants.hpp"

namespace SolTrace::Data {

int st_sun_position(double lat, double day, double hour,
                    double *x, double *y, double *z)
{
    /*
    computes the sun vector xyz given arguments
    lat : [deg] latitude
    day : [] day of the year
    hour : [hour] solar time. 12.00 corresponds to sun at maximum elevation and does not necessarily match local time

    xyz coordinate system:
        x: +west
        y: +zenith
        z: +north
    */
    double Elevation, Azimuth, Zenith;

    legacy_sun_position(lat, day, hour, &Azimuth, &Zenith);
    Elevation = 90.0 - Zenith;

    // TODO: Update coordinate system
    *x = -sin(Azimuth * D2R) * cos(Elevation * D2R);
    *y = sin(Elevation * D2R);
    *z = cos(Azimuth * D2R) * cos(Elevation * D2R);

    return 1;
}

void legacy_sun_position(double lat, double day, double hour, double* azimuth, double* zenith)
{
    /*
    Computes the sun Azimuth and Zenith angles using the legacy SolTrace method. This method is similar to that presented in
    Duffie, J.A. and Beckman, W.A., "Solar Engineering of Thermal Processes," 4th Edition, Wiley, 2013. but not identical.
        lat : [deg] latitude
        day : [] day of the year
        hour : [hour] solar time. 12.00 corresponds to sun at maximum elevation and does not necessarily match local time

    Returns:
        azimuth : [deg] azimuth angle, measured east from north
        zenith : [deg] zenith angle

    */

    double Declination, HourAngle, Elevation, Azimuth;

    Declination = R2D * asin(0.39795 * cos(0.98563 * D2R * (day - 173)));
    HourAngle = 15 * (hour - 12);
    Elevation = R2D * asin(sin(Declination * D2R) * sin(lat * D2R) + cos(Declination * D2R) * cos(HourAngle * D2R) * cos(lat * D2R));
    Azimuth = R2D * acos((sin(D2R * Declination) * cos(D2R * lat) - cos(D2R * Declination) * sin(D2R * lat) * cos(D2R * HourAngle)) / cos(D2R * Elevation) + 0.0000000001);
    if (sin(HourAngle * D2R) > 0.0)
        Azimuth = 360 - Azimuth;

    *azimuth = Azimuth;
    *zenith = 90.0 - Elevation;
}

void duffie_sun_position(double lat, double lng, double tz, double day, double hour, double* azimuth, double* zenith)
{
    /*
    Computes the sun Azimuth and Zenith angles using the method presented in
    Duffie, J.A. and Beckman, W.A., "Solar Engineering of Thermal Processes," 4th Edition, Wiley, 2013.

    lat : [deg] latitude
    day : [] day of the year
    hour : [hour] solar time. 12.00 corresponds to sun at maximum elevation and does not necessarily match local time

    Returns:
    azimuth : [deg] azimuth angle, measured east from north
    zenith : [deg] zenith angle
    */
    // TODO: we could improve this by 
    //  - calculate a fractional day of year
    //  - Use More accurate equation for declination


    double B = (day - 1) * 360. / 365.; // [degrees] Equation 1.4.2
    double E = 229.2 * (0.000075 + 0.001868 * cos(D2R * B) - 0.032077 * sin(D2R * B) - 0.014615 * cos(D2R * 2 * B) - 0.04089 * sin(D2R * 2 * B)); // [minutes] Equation 1.5.3 
    double solar_time = hour + (lng / 15.0 - tz) + E / 60.0;  // [hour] solar time, Equation 1.5.2

    double HourAngle = 15 * (solar_time - 12);
    //double Declination = 23.45 * sin(D2R * (360.0 * (284 + day) / 365.0));
    double Declination = R2D * (0.006918 - 0.399912 * cos(D2R * B) + 0.070257 * sin(D2R * B) - 0.006758 * cos(D2R * 2 * B) + 0.000907 * sin(D2R * 2 * B) - 0.002697 * cos(D2R * 3 * B) + 0.00148 * sin(D2R * 3 * B));
    double Zenith = R2D * acos(cos(D2R * lat) * cos(D2R * Declination) * cos(D2R * HourAngle) + sin(D2R * lat) * sin(D2R * Declination));
    double ratio = (cos(D2R * Zenith) * sin(D2R * lat) - sin(D2R * Declination)) / (sin(D2R * Zenith) * cos(D2R * lat));

    double Azimuth;
    // adjusting numerical error
    if (ratio > 1.0 && (ratio - 1.0) < 1.e-6)
        Azimuth = 0.0;
    else if (ratio < -1.0 && (-1.0 - ratio) < 1.e-6)
        Azimuth = 180.0;
    else 
        Azimuth = R2D * fabs(acos(ratio));

    if (HourAngle < 0)
        Azimuth *= -1.0;
    Azimuth += 180.0;    // Shift to typical North reference (east positive angle)

    *azimuth = Azimuth;
    *zenith = Zenith;
}



} // namespace SolTrace::Data

1
2	#include "basic_sun_position.hpp"
3	#include "constants.hpp"
4
5	namespace SolTrace::Data {
6
7	int st_sun_position(double lat, double day, double hour,	2✔
8	double x, double y, double *z)
9	{
10	/*
11	computes the sun vector xyz given arguments
12	lat : [deg] latitude
13	day : [] day of the year
14	hour : [hour] solar time. 12.00 corresponds to sun at maximum elevation and does not necessarily match local time
15
16	xyz coordinate system:
17	x: +west
18	y: +zenith
19	z: +north
20	*/
21	double Elevation, Azimuth, Zenith;
22
23	legacy_sun_position(lat, day, hour, &Azimuth, &Zenith);	2✔
24	Elevation = 90.0 - Zenith;	2✔
25
26	// TODO: Update coordinate system
27	x = -sin(Azimuth D2R) * cos(Elevation * D2R);	2✔
28	y = sin(Elevation D2R);	2✔
29	z = cos(Azimuth D2R) * cos(Elevation * D2R);	2✔
30
31	return 1;	2✔
32	}
33
34	void legacy_sun_position(double lat, double day, double hour, double* azimuth, double* zenith)	4✔
35	{
36	/*
37	Computes the sun Azimuth and Zenith angles using the legacy SolTrace method. This method is similar to that presented in
38	Duffie, J.A. and Beckman, W.A., "Solar Engineering of Thermal Processes," 4th Edition, Wiley, 2013. but not identical.
39	lat : [deg] latitude
40	day : [] day of the year
41	hour : [hour] solar time. 12.00 corresponds to sun at maximum elevation and does not necessarily match local time
42
43	Returns:
44	azimuth : [deg] azimuth angle, measured east from north
45	zenith : [deg] zenith angle
46
47	*/
48
49	double Declination, HourAngle, Elevation, Azimuth;
50
51	Declination = R2D * asin(0.39795 * cos(0.98563 * D2R * (day - 173)));	4✔
52	HourAngle = 15 * (hour - 12);	4✔
53	Elevation = R2D * asin(sin(Declination * D2R) * sin(lat * D2R) + cos(Declination * D2R) * cos(HourAngle * D2R) * cos(lat * D2R));	4✔
54	Azimuth = R2D * acos((sin(D2R * Declination) * cos(D2R * lat) - cos(D2R * Declination) * sin(D2R * lat) * cos(D2R * HourAngle)) / cos(D2R * Elevation) + 0.0000000001);	4✔
55	if (sin(HourAngle * D2R) > 0.0)	4✔
56	Azimuth = 360 - Azimuth;	1✔
57
58	*azimuth = Azimuth;	4✔
59	*zenith = 90.0 - Elevation;	4✔
60	}	4✔
61
62	void duffie_sun_position(double lat, double lng, double tz, double day, double hour, double* azimuth, double* zenith)	541✔
63	{
64	/*
65	Computes the sun Azimuth and Zenith angles using the method presented in
66	Duffie, J.A. and Beckman, W.A., "Solar Engineering of Thermal Processes," 4th Edition, Wiley, 2013.
67
68	lat : [deg] latitude
69	day : [] day of the year
70	hour : [hour] solar time. 12.00 corresponds to sun at maximum elevation and does not necessarily match local time
71
72	Returns:
73	azimuth : [deg] azimuth angle, measured east from north
74	zenith : [deg] zenith angle
75	*/
76	// TODO: we could improve this by
77	// - calculate a fractional day of year
78	// - Use More accurate equation for declination
79
80
81	double B = (day - 1) * 360. / 365.; // [degrees] Equation 1.4.2	541✔
82	double E = 229.2 * (0.000075 + 0.001868 * cos(D2R * B) - 0.032077 * sin(D2R * B) - 0.014615 * cos(D2R * 2 * B) - 0.04089 * sin(D2R * 2 * B)); // [minutes] Equation 1.5.3	541✔
83	double solar_time = hour + (lng / 15.0 - tz) + E / 60.0; // [hour] solar time, Equation 1.5.2	541✔
84
85	double HourAngle = 15 * (solar_time - 12);	541✔
86	//double Declination = 23.45 * sin(D2R * (360.0 * (284 + day) / 365.0));
87	double Declination = R2D * (0.006918 - 0.399912 * cos(D2R * B) + 0.070257 * sin(D2R * B) - 0.006758 * cos(D2R * 2 * B) + 0.000907 * sin(D2R * 2 * B) - 0.002697 * cos(D2R * 3 * B) + 0.00148 * sin(D2R * 3 * B));	541✔
88	double Zenith = R2D * acos(cos(D2R * lat) * cos(D2R * Declination) * cos(D2R * HourAngle) + sin(D2R * lat) * sin(D2R * Declination));	541✔
89	double ratio = (cos(D2R * Zenith) * sin(D2R * lat) - sin(D2R * Declination)) / (sin(D2R * Zenith) * cos(D2R * lat));	541✔
90
91	double Azimuth;
92	// adjusting numerical error
93	if (ratio > 1.0 && (ratio - 1.0) < 1.e-6)	541✔
NEW 94	Azimuth = 0.0;	×
95	else if (ratio < -1.0 && (-1.0 - ratio) < 1.e-6)	541✔
NEW 96	Azimuth = 180.0;	×
97	else
98	Azimuth = R2D * fabs(acos(ratio));	541✔
99
100	if (HourAngle < 0)	541✔
101	Azimuth *= -1.0;	239✔
102	Azimuth += 180.0; // Shift to typical North reference (east positive angle)	541✔
103
104	*azimuth = Azimuth;	541✔
105	*zenith = Zenith;	541✔
106	}	541✔
107
108
109
110	} // namespace SolTrace::Data

NREL / SolTrace / 19315321183

Source File Press 'n' to go to next uncovered line, 'b' for previous

Source File
Press 'n' to go to next uncovered line, 'b' for previous